Smoking and COVID Living Review

0.1 Data preparation

0.1.1 Load packages

• Load packages: load packages and useful scripts for the analysis

library('tidyverse', quietly = T)
library('ggplot2', quietly = T)
library('dplyr', quietly = T)
library('readxl', quietly = T)
library('linelist', quietly = T)
library('googlesheets4', quietly = T)
library('naniar', quietly = T)
library('gridExtra', quietly = T)
library('snakecase', quietly = T)
library('epitools', quietly = T)
library('meta', quietly = T)
library('knitr', quietly = T)
library('kableExtra', quietly = T)
library('snakecase', quietly = T)
library('here', quietly = T)
library('flextable', quietly = T)

0.1.2 Download from Google sheets

• Load data: import raw dataset directly from Google sheets

here()

## [1] "C:/Users/ucbtds4/Google Drive/PhD/COVID-19/Smoking and COVID data/smoking_covid/smoking_covid"

sheets_id <- as_sheets_id('https://docs.google.com/spreadsheets/d/15avypGR8ypJngWQEmFIzFrOOwXPY3xezUHQU6jgV7d0/edit?usp=sharing')
search_details <- read_sheet(sheets_id, range = 'article_screening')
data_study_general <- read_sheet(sheets_id, range =  'general_details')
table_1 <- read_sheet(sheets_id, sheet = 'pop_descriptives')
table_2 <-  read_sheet(sheets_id, sheet = 'testing')
table_3 <-  read_sheet(sheets_id, sheet = 'hospitalisation')
table_4 <-  read_sheet(sheets_id, sheet = 'severity')
table_5 <-  read_sheet(sheets_id, sheet = 'mortality')
table_6 <-  read_sheet(sheets_id, sheet = 'quality_appraisal')

0.1.3 Clean data and save updated versions as RDS files

search_details <- search_details %>%
  clean_data() %>%
  rename(date_screened = date_screening_in_format_2020_05_22) %>%
  write_rds(here::here('data_clean', 'search_details.rds'))


protect_columns_1 <- names(data_study_general) %in% 'doi'

data_study_general <- data_study_general %>%
  rename('notes' = 11, date_published = `date_published_format(2020-05-12)`) %>%
  clean_data(protect = protect_columns_1) %>%
  mutate(study_id = 1:length(lead_author)) %>%
  write_rds(here::here('data_clean', 'data_study_general.rds'))

review_details <- data_study_general %>%
  select(lead_author, date_published, country, review_version, study_id) %>%
  write_rds(here::here('data_clean', 'review_details.rds'))

table_1 <- table_1 %>%
  clean_data() %>%
  mutate(., lower_range = sub('\\_.*', '', .$range)) %>%
  mutate(., upper_range = sub('.*_', '',.$range,)) %>%
  select(1:9, 22:23, 11:21) %>%
  filter(!is.na(lead_author)) %>%
  left_join(., review_details, by = 'lead_author') %>%
  write_rds(here::here('data_clean', 'table_1.rds'))

table_2 <-  table_2 %>%
  clean_data() %>%
  filter(data_on_testing == TRUE) %>%
  left_join(., review_details, by = 'lead_author') %>%
  write_rds(here::here('data_clean', 'table_2.rds'))

table_3 <-  table_3 %>%
  clean_data() %>%
  filter(data_on_hospitalisation == TRUE) %>%
  left_join(., review_details, by = 'lead_author') %>%
  write_rds(here::here('data_clean', 'table_3.rds'))

table_4 <-  table_4 %>%
  clean_data() %>%
  filter(data_disease_severity == TRUE) %>%
  left_join(., review_details, by = 'lead_author') %>%
  write_rds(here::here('data_clean', 'table_4.rds'))

table_5 <- table_5 %>%
  clean_data() %>%
  filter(data_on_deaths == TRUE) %>%
  left_join(., review_details, by = 'lead_author') %>%
  write_rds(here::here('data_clean', 'table_5.rds'))


protect_columns_2 <- !names(table_6) %in% 'Author'

table_6 <-  table_6 %>%
  clean_data(protect = protect_columns_2) %>%
  left_join(., table_1 %>%
              select(lead_author,
                     not_stated,
                     missing,
                     total) %>%
              rename('author' = lead_author),
            by = 'author') %>%
  mutate(missingness = rowSums(.[17:18], na.rm = T),
         missingness_percentage = (missingness/total)*100) %>%
  select(-c(not_stated, missing, total, missingness)) %>%
  write_rds(here::here('data_clean', 'table_6.rds'))

a <- data_study_general %>%
  select(lead_author, date_published, source, study_id) %>%
  rename('Lead Author' = lead_author,
         'Date Published' = date_published,
         'Publication Source' = source,
         'Study ID' = study_id)

a$`Lead Author` <- to_upper_camel_case(a$`Lead Author`, sep_out = ", ")
a$`Publication Source` <- to_title_case(a$`Publication Source`)
a$`Publication Source` <- if_else(str_length(a$`Publication Source`) < 5,
                                  toupper(a$`Publication Source`),
                                  to_title_case(a$`Publication Source`))

a <- flextable(a)%>%
  set_caption(a, caption = 'All studies included in data extraction') %>%
  set_table_properties(width = 1, layout = 'autofit') %>%
  save_as_docx(a, path = here('data_clean', 'all_studies.docx'))

All data has now been stored in .rds format in the /data_clean folder. This can be analysed outside of the report using pre-generated scripts or user created scripts.

0.2 Preparing the descriptive analysis

0.2.1 Choosing which data we want to include in the current review?

date_of_update <- Sys.Date()

prev_versions <- c('v1', 'v2', 'v3')
#This can categorise which studies we want to look at
current_version <- c('v4')
#This will categorise which studies we are including in the current report

analysed_versions <- c('v1', 'v2', 'v3', 'v4')
#This will incorporate all studies into the current version of the report

exclude_from_analysis <- c('isaric_1', 'isaric_2', 'isaric_3', 'miyara_old', 'isaric_4', 'mehra')

exclude_from_qa <- c('isaric_1', 'isaric_2', 'isaric_3', 'miyara', 'isaric_4', 'mehra')

We are currently on Version v4. We have excluded from the analysis the following studies: isaric_1, isaric_2, isaric_3, miyara_old, isaric_4, and mehra

0.2.2 Filter the included studies

data_study_general <- data_study_general %>%
  filter(review_version %in% analysed_versions) %>%
  filter(!(lead_author %in% exclude_from_analysis)) %>%
  write_csv(., here::here('data_clean', 'data_study_general.csv'))

table_1 <- table_1 %>%
  filter(review_version %in% analysed_versions) %>%
  select(-review_version) %>%
  filter(!(lead_author %in% exclude_from_analysis))

table_2 <- table_2 %>%
  filter(review_version %in% analysed_versions) %>%
  select(-review_version) %>%
  filter(!(lead_author %in% exclude_from_analysis))

table_3 <- table_3 %>%
  filter(review_version %in% analysed_versions) %>%
  select(-review_version) %>%
  filter(!(lead_author %in% exclude_from_analysis))

table_4 <- table_4 %>%
  filter(review_version %in% analysed_versions) %>%
  select(-review_version) %>%
  filter(!(lead_author %in% exclude_from_analysis))

table_5 <- table_5 %>%
  filter(review_version %in% analysed_versions) %>%
  select(-review_version) %>%
  filter(!(lead_author %in% exclude_from_analysis))

table_6 <- table_6 %>%
  filter(!(author %in% exclude_from_qa))

a <- data_study_general %>%
  select(lead_author, date_published, source, study_id) %>%
  rename('Lead Author' = lead_author,
         'Date Published' = date_published,
         'Publication Source' = source,
         'Study ID' = study_id)

a$`Lead Author` <- to_upper_camel_case(a$`Lead Author`, sep_out = ", ")
a$`Publication Source` <- to_title_case(a$`Publication Source`)
a$`Publication Source` <- if_else(str_length(a$`Publication Source`) < 5,
                                  toupper(a$`Publication Source`),
                                  to_title_case(a$`Publication Source`))

a <- flextable(a)%>%
  set_table_properties(width = 1, layout = 'autofit') %>%
  set_caption(caption = 'Studies included in the current analysis')

save_as_docx(a, path = here('data_clean', 'included_studies.docx'))

a

Studies included in the current analysis

Lead Author	Date Published	Publication Source	Study ID
Guan, Ni	2020-02-28	NEJM	1
Guan, Liang	2020-03-26	European Respiratory Journal	2
Lian	2020-03-25	Clinical Infectious Diseases	3
Jin	2020-03-24	Gastroduodenal	4
Chen	2020-03-26	BMJ	5
Zhou, Yu	2020-03-11	The Lancet	6
Mo	2020-03-16	Clinical Infectious Diseases	7
Zhang, Dong	2020-02-19	Allergy	8
Wan	2020-03-21	Journal of Medical Virology	9
Liu, Tao	2020-02-28	Chinese Medical Journal	10
Huang, Wang	2020-01-24	The Lancet	11
Zhang, Cai	2020-03-20	International Journal of Infectious Diseases	12
Guo	2020-03-27	Jama Cardiology	13
Liu, Ming	2020-03-12	Medrxiv	14
Huang, Yang	2020-03-05	Medrxiv	15
Xu	2020-03-08	Medrxiv	16
Li	2020-02-12	Medrxiv	17
Rentsch	2020-04-14	Medrxiv	18
Hu	2020-03-25	Medrxiv	19
Wang, Pan	2020-03-24	International Journal of Infectious Diseases	20
Petrilli	2020-04-11	Medrxiv	22
Chow, Us, Cdc	2020-03-31	Report	23
Dong, Cao	2020-03-20	Allergy	25
Kim	2020-04-01	Journal of Korean Medical Science	26
Shi, Yu	2020-03-18	Bmc Critical Care	27
Yang, Yu	2020-02-24	The Lancet Respiratory Medicine	28
Argenziano	2020-04-22	Medrxiv	29
Solis	2020-04-25	Medrxiv	30
Richardson	2020-04-22	JAMA	31
Fontanet	2020-04-23	Medrxiv	32
Zheng, Gao	2020-04-19	Metabolism	33
Liao, Feng	2020-04-24	Medrxiv	34
Rodriguez	2020-04-24	Medrxiv	35
Magagnoli	2020-04-16	Medrxiv	36
Shi, Ren	2020-04-23	Medrxiv	37
Hadjadj	2020-04-23	Medrxiv	38
Niedzwiedz	2020-04-30	Medrxiv	39
Gold, Us, Cdc	2020-04-20	Cdc Report	41
Yu, Cai	2020-04-27	Clinical Therapeutics	44
Zheng, Xiong	2020-04-30	Pharmacological Research	45
Miyara, Updated	2020-05-09	Qeios	46
De, La, Rica	2020-05-11	Medrxiv	47
Yin, Yang	2020-05-10	Medrxiv	48
Gaibazzi	2020-05-10	Medrxiv	49
Shi, Zuo	2020-05-10	Medrxiv	50
Cho	2020-05-11	Medrxiv	51
Allenbach	2020-05-08	Medrxiv	52
Robilotti	2020-05-08	Medrxiv	53
Opensafely, Collaborative	2020-05-07	Medrxiv	54
Borobia	2020-05-06	Medrxiv	55
Giacomelli	2020-05-06	Medrxiv	56
Shah	2020-05-06	Medrxiv	57
Bello, Chavolla	2020-05-06	Medrxiv	58
Kolin	2020-05-05	Medrxiv	59
Lubetzky	2020-05-08	Medrxiv	60
Goyal	2020-04-17	NEJM	61
Feng	2020-04-10	American Journal of Respiratory and Critical Care Medicine	62
Yao	2020-04-24	Polish Archives of Internal Medicine	63
Sami	2020-05-15	Medrxiv	64
Almazeedi	2020-05-15	Medrxiv	65
Carrillo, Vega	2020-05-14	Medrxiv	66
Yanover	2020-05-13	Medrxiv	67
Hamer	2020-05-13	Medrxiv	68
Regina	2020-05-14	Medrxiv	69
De, Lusignan	2020-05-15	Lancet Infectious Diseases	71
Targher	2020-05-13	Diabetes Metabolism	72
Valenti	2020-05-18	Medrxiv	73
Feuth	2020-05-18	Medrxiv	74
Heng	2020-05-18	Medrxiv	75
Parrotta	2020-05-18	Medrxiv	76
Shekar	2020-05-18	Medrxiv	77
Mejia, Vilet	2020-05-16	Medrxiv	78
Chen, Jiang	2020-05-16	Medrxiv	79
Jiong	2020-05-16	Medrxiv	80
Rimland	2020-05-19	Medrxiv	81
Palaiodimos	2020-05-15	Metabolism	82
Ip	2020-05-25	Medrxiv	83
Heili, Frades	2020-05-25	Medrxiv	84
Vaquero, Roncero	2020-05-24	Medrxiv	85
Kim, Garg	2020-05-22	Medrxiv	86
Wu	2020-05-21	Medrxiv	87
Hopkinson	2020-05-21	Medrxiv	88
Shi, Zhao	2020-05-20	Medrxiv	89
Kimmig	2020-05-20	Medrxiv	90
Al, Hindawi	2020-05-20	Medrxiv	91
Basse	2020-05-19	Medrxiv	92
Freites	2020-05-19	Medrxiv	93
Alshami	2020-05-19	Medrxiv	94
Russell	2020-05-19	Medrxiv	95
Berumen	2020-05-26	Medrxiv	96
Gianfrancesco	2020-05-29	Annals of the Rheumatic Diseases	97
Li, Long	2020-05-28	Movement Disorders	98
Batty	2020-06-01	Medrxiv	99
Israel	2020-06-01	Medrxiv	100
Del, Valle	2020-05-30	Medrxiv	101
Zuo, Zuo	2020-05-29	Medrxiv	102
Chaudhry	2020-05-29	Medrxiv	103
Louis	2020-06-28	Medrxiv	104
Soto, Mota	2020-05-27	Medrxiv	105
Patel	2020-05-26	Medrxiv	106
Garibaldi	2020-05-26	Medrxiv	107
Docherty	2020-05-22	BMJ	108
Boulware	2020-03-06	NEJM	109

0.3 Descriptive analysis

0.3.1 The PRISMA diagram will be updated with these numbers

#Number screened to add to PRISMA
from_prev_version <- data_study_general %>%
  filter(review_version %in% prev_versions) %>%
  filter(!(lead_author %in% exclude_from_analysis)) %>%
  tally()

source(here::here('scripts', 'prisma_function.R'))

PRISMA_v4 <- PRISMA(search_details)
PRISMA_v4 <- as.data.frame(PRISMA_v4)
a <- PRISMA_v4 %>%
  t() %>%
  as.data.frame() %>%
  add_rownames() %>%
  flextable() %>%
  set_header_labels(rowname = 'PRISMA category',
                    V1 = 'Number of studies') %>%
  set_table_properties(width = 1, layout = 'autofit') %>%
  set_caption(caption = 'Data for PRISMA')

save_as_docx(a, path = here('data_clean', 'PRISMA_details.docx'))
a

Data for PRISMA

PRISMA category	Number of studies
ovid_number_results	297
medrxiv_number_results	158
other_source_included	2
titles_screened	457
excluded_abstract	249
full_texts_reviewed	208
full_texts_excluded	170
full_texts_included_narrative	38
total_all_versions	103

0.3.2 We summarise where the studies were conducted

#Countries
country <- table_1$country %>%
  to_upper_camel_case() %>%
  recode(., 'Usa' = 'USA', 'Uk' = 'UK', 'SaudiArabia' = 'Saudi Arabia') %>%
  table() %>%
  as.data.frame() %>%
  arrange(desc(Freq)) %>%
  rename('Country' = 1,
         'Number' = 2)

ggplot(country, aes(x = reorder(Country, desc(-Number)), y = Number))+
  geom_col()+
  coord_flip()+
  theme_bw()+
  labs(title = 'Countries where studies were performed',
       y = 'Number of studies',
       x = 'Country')

a <- flextable(country) %>%
  set_table_properties(width = 0.5, layout = 'autofit')
a

Country	Number
China	36
USA	25
UK	10
Mexico	6
Spain	6
France	5
Italy	4
Multiple	3
Israel	2
Finland	1
Iran	1
Korea	1
Kuwait	1
Saudi Arabia	1
Switzerland	1

save_as_docx(a, path = here('data_clean', 'Countries_performing_studies.docx'))  

0.3.3 We sumarise the settings in which they were conducted

#Setting
setting <- to_upper_camel_case(data_study_general$study_setting, sep_out = ' ') %>%
  table(.) %>%
  sort('Number of studies', decreasing = T) %>%
  as.data.frame() %>%
  arrange(desc(Freq)) %>%
  rename('Setting' = 1,
         'Number' = 2)

a <- flextable(setting) %>%
  set_table_properties(width = 0.5, layout = 'autofit')
a

Setting	Number
Hospital	81
Community And Hospital	16
Community	4
Not Stated	1
Quarantine Centre	1

save_as_docx(a, path = here('data_clean', 'Countries_performing_studies.docx'))  

0.3.4 We update the number of participants we have included

• Summary statistics: This table will give us the mean, median and IQR of the samples in the review

#Numbers
summary(table_1$sample_size)

##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
##        7       96      316   199095     1052 17425445

• The total number of participants

#Number of participants
sum(table_1$sample_size)

## [1] 20506745

0.3.5 We want to know if they use a CRF or electronic health record

#Source of data
table_1$data_source %>%
  recode('electronic_health_records' = 'Electronic health record',
         'case_report_form' = 'Case report form',
         'not_stated' = 'Not stated') %>%
  table() %>%
  as.data.frame() %>%
  arrange(desc(Freq)) %>%
  rename('Data source' = 1,
         'Number' = 2) %>%
  flextable() %>%
  set_table_properties(width = 0.5, layout = 'autofit')

Data source	Number
Electronic health record	62
Case report form	26
Not stated	15

0.3.6 Now about the smoking

• Smokers

#Studies reporting smokers
current_smok <- table_1 %>%
  filter(., current_smoker != 'NA')
nrow(current_smok)

## [1] 69

• Former smokers

#Studies reporting former smokers
former_smok <- table_1 %>%
  filter(., former_smoker != 'NA')
nrow(former_smok)

## [1] 31

• Never smokers

#Studies reporting never smokers
never_smok <- table_1 %>%
  filter(., never_smoker != 'NA')
nrow(never_smok)

## [1] 29

• Current/former smokers

#Studies reporting current/former smokers
current_former_smok <- table_1 %>%
  filter(., current_former_smoker != 'NA')
nrow(current_former_smok)

## [1] 34

• Missing, never/unknown and not stated

#Studies reporting missing data
missing_smok <- table_1 %>%
  filter(., missing != 'NA')
nrow(missing_smok)

## [1] 27

#Studies reporting never/unknown
never_smok_unknown <- table_1 %>%
  filter(., never_smoker_unknown != 'NA')
nrow(never_smok_unknown)

## [1] 8

#Studies with not stated
smok_not_stated <- table_1 %>%
  filter(., not_stated != 'NA')
nrow(smok_not_stated)

## [1] 65

• Studies with complete smoking status

#Studies reporting current, former and never smoking status
full_smoking_status <- table_1 %>%
  filter(lead_author %in% current_smok$lead_author) %>%
  filter(lead_author %in% former_smok$lead_author) %>%
  filter(lead_author %in% never_smok$lead_author)

b <- full_smoking_status %>%
  select(lead_author, sample_size, current_smoker, former_smoker, never_smoker, study_id) %>%
  rename('Lead author' = lead_author,
         'Sample size' = sample_size,
         'Current smokers' = current_smoker,
         'Former smokers' = former_smoker,
         'Never smokers' = never_smoker,
         'Study ID' = study_id)

b$`Lead author` <- to_upper_camel_case(b$`Lead author`, sep_out = ", ")

b <- flextable(b) %>%
  set_table_properties(width = 1, layout = 'autofit') %>%
  set_caption(caption = 'Studies with complete smoking status data')
b

Studies with complete smoking status data

Lead author	Sample size	Current smokers	Former smokers	Never smokers	Study ID
Guan, Ni	1099	137	21	927	1
Rentsch	3528	159	179	216	18
Argenziano	1000	49	179	772	29
Hadjadj	50	1	9	40	38
Niedzwiedz	1474	42745	147235	235711	39
Miyara, Updated	479	32	152	286	46
Gaibazzi	441	21	44	376	49
Cho	1331	253	359	719	51
Robilotti	423	9	159	248	53
Opensafely, Collaborative	17425445	2962373	5737545	8000204	54
Shah	316	52	56	133	57
Kolin	1474	213	592	657	59
Yanover	4353	514	129	3710	67
Hamer	387109	37426	134855	214828	68
De, Lusignan	3802	413	1753	1125	71
Feuth	28	3	8	17	74
Parrotta	76	2	20	52	76
Ip	2512	77	448	1620	83
Al, Hindawi	31	1	22	8	91
Russell	106	7	26	42	95
Israel	24087	4114	3141	16832	100
Docherty	20133	852	4364	8968	108

nrow(full_smoking_status)

## [1] 22

#These are the studies with data on current, former and never smokers

#Studies reporting current or current/former and never smoking
semi_full_smoking_status <- table_1 %>%
  filter(lead_author %in% current_former_smok$lead_author) %>%
  filter(lead_author %in% never_smok$lead_author) %>%
  filter(!lead_author %in% full_smoking_status$lead_author)

b <- semi_full_smoking_status %>%
  select(lead_author, sample_size, current_former_smoker, never_smoker, study_id) %>%
  rename('Lead author' = lead_author,
         'Sample size' = sample_size,
         'Current/former smokers' = current_former_smoker,
         'Never smokers' = never_smoker,
         'Study ID' = study_id)

b$`Lead author` <- to_upper_camel_case(b$`Lead author`, sep_out = ", ")

b <- flextable(b) %>%
  set_table_properties(width = 1, layout = 'autofit') %>%
  set_caption(caption = 'Studies with partially complete smoking status (current/former rather than current and former)')
b

Studies with partially complete smoking status (current/former rather than current and former)

Lead author	Sample size	Current/former smokers	Never smokers	Study ID
Guan, Liang	1590	111	1479	2
Richardson	5700	558	3009	31
Rodriguez	7	3	4	35
Zheng, Xiong	73	8	65	45
Giacomelli	233	70	163	56
Palaiodimos	200	65	135	82
Gianfrancesco	600	129	389	97

nrow(semi_full_smoking_status)

## [1] 7

#Remaining studies
incomplete_smoking_status <- table_1 %>%
  filter(!lead_author %in% full_smoking_status$lead_author) %>%
  filter(!lead_author %in% semi_full_smoking_status$lead_author)

b <- incomplete_smoking_status %>%
  select(lead_author, sample_size, current_smoker, current_former_smoker, former_smoker, never_smoker, not_stated, missing, study_id) %>%
  rename('Lead author' = lead_author,
         'Sample size' = sample_size,
         'Current smoker' = current_smoker,
         'Current/former smokers' = current_former_smoker,
         'Former smoker' = former_smoker,
         'Never smokers' = never_smoker,
         'Not stated' = not_stated,
         'Missing' = missing,
         'Study ID' = study_id)

b$`Lead author` <- to_upper_camel_case(b$`Lead author`, sep_out = ", ")

b <- flextable(b) %>%
  set_table_properties(width = 1, layout = 'autofit') %>%
  set_caption(caption = 'Studies with incomplete smoking status')
b

Studies with incomplete smoking status

Lead author	Sample size	Current smoker	Current/former smokers	Former smoker	Never smokers	Not stated	Missing	Study ID
Lian	788	54	NA	NA	NA	734	NA	3
Jin	651	41	NA	NA	NA	610	NA	4
Chen	548	24	NA	14	NA	510	NA	5
Zhou, Yu	191	11	NA	NA	NA	180	NA	6
Mo	155	6	NA	NA	NA	149	NA	7
Zhang, Dong	140	2	NA	7	NA	131	NA	8
Wan	135	9	NA	NA	NA	126	NA	9
Liu, Tao	78	NA	5	NA	NA	73	NA	10
Huang, Wang	41	3	NA	NA	NA	38	NA	11
Zhang, Cai	645	41	NA	NA	NA	604	NA	12
Guo	187	18	NA	NA	NA	169	NA	13
Liu, Ming	41	4	NA	NA	NA	37	NA	14
Huang, Yang	36	NA	4	NA	NA	32	NA	15
Xu	53	6	NA	NA	NA	47	NA	16
Li	17	3	NA	NA	NA	14	NA	17
Hu	323	NA	38	NA	NA	285	NA	19
Wang, Pan	125	NA	16	NA	NA	109	NA	20
Petrilli	4103	212	NA	666	NA	NA	NA	22
Chow, Us, Cdc	7162	96	NA	165	NA	NA	6901	23
Dong, Cao	9	1	NA	NA	NA	8	NA	25
Kim	28	5	NA	NA	NA	22	1	26
Shi, Yu	487	NA	40	NA	NA	434	13	27
Yang, Yu	52	2	NA	NA	NA	50	NA	28
Solis	650	61	NA	NA	NA	589	NA	30
Fontanet	661	69	NA	NA	NA	NA	NA	32
Zheng, Gao	66	8	NA	NA	NA	58	NA	33
Liao, Feng	1848	NA	8	NA	NA	140	1700	34
Magagnoli	368	NA	52	NA	NA	316	NA	36
Shi, Ren	134	NA	14	NA	NA	120	NA	37
Gold, Us, Cdc	305	16	NA	NA	NA	NA	289	41
Yu, Cai	95	8	NA	NA	NA	87	NA	44
De, La, Rica	48	NA	10	NA	NA	37	1	47
Yin, Yang	106	NA	18	NA	NA	88	NA	48
Shi, Zuo	96	NA	29	NA	NA	67	NA	50
Allenbach	152	NA	10	NA	NA	142	NA	52
Borobia	2226	157	NA	NA	NA	2069	NA	55
Bello, Chavolla	62489	NA	6209	NA	NA	56280	NA	58
Lubetzky	54	NA	12	NA	NA	42	NA	60
Goyal	393	20	NA	NA	NA	373	NA	61
Feng	476	44	NA	NA	NA	410	22	62
Yao	108	4	NA	NA	NA	NA	104	63
Sami	490	69	NA	NA	NA	NA	NA	64
Almazeedi	1096	44	NA	NA	NA	NA	NA	65
Carrillo, Vega	10544	936	NA	NA	NA	9608	NA	66
Regina	200	9	NA	NA	NA	191	NA	69
Targher	339	28	NA	NA	NA	311	NA	72
Valenti	789	204	NA	NA	NA	585	NA	73
Heng	51	7	NA	NA	NA	44	NA	75
Shekar	50	24	NA	NA	NA	26	NA	77
Mejia, Vilet	329	NA	23	NA	NA	306	NA	78
Chen, Jiang	135	NA	13	NA	NA	122	NA	79
Jiong	1008	57	NA	NA	NA	484	467	80
Rimland	11	1	NA	NA	NA	9	NA	81
Heili, Frades	4712	233	NA	820	NA	NA	526	84
Vaquero, Roncero	146	NA	10	NA	NA	136	NA	85
Kim, Garg	2491	150	NA	642	NA	NA	2	86
Wu	174	NA	58	NA	NA	116	NA	87
Hopkinson	2401982	NA	263257	NA	NA	2138725	NA	88
Shi, Zhao	101	NA	5	NA	NA	96	NA	89
Kimmig	60	NA	22	NA	NA	38	NA	90
Basse	141	25	NA	NA	NA	116	NA	92
Freites	123	4	NA	NA	NA	119	NA	93
Alshami	128	20	NA	3	NA	105	NA	94
Berumen	102875	NA	9921	NA	NA	NA	NA	96
Li, Long	145	NA	8	NA	NA	137	NA	98
Batty	908	102	NA	NA	NA	806	NA	99
Del, Valle	1484	82	NA	346	NA	680	376	101
Zuo, Zuo	44	NA	12	NA	NA	32	NA	102
Chaudhry	40	NA	6	NA	NA	33	1	103
Louis	22	NA	10	NA	NA	12	NA	104
Soto, Mota	400	NA	48	NA	NA	352	NA	105
Patel	104	43	NA	NA	NA	NA	10	106
Garibaldi	832	46	NA	188	NA	598	NA	107
Boulware	821	27	NA	NA	NA	794	NA	109

nrow(incomplete_smoking_status)

## [1] 74

0.4 Smoking prevalence by country

#Smoking prevalence by country
country_prevalence_list <- table_1 %>%
  group_by(country) %>%
  mutate(., current_smok_percentage = current_smoker/total*100) %>%
  mutate(., former_smok_percentage = former_smoker/total*100) %>%
  mutate(., missing_percentage = missing/total*100) %>%
  mutate(., not_stated_percentage = not_stated/total*100) %>%
  mutate(never_smoker_percentage = never_smoker/total*100) %>%
  group_split(country)

We have split the data into a list where each country can be analysed separately if required

0.4.1 Graph for smoking prevalence in the countries included

Country smoking prevalence

0.5 Final results tables and meta-analyses

0.5.1 Table 1: Included studies

#Updating table 1
table_1_word <- table_1 %>%
  mutate(., current_percentage = current_smoker/total*100) %>%
  mutate(., former_percentage = former_smoker/total*100) %>%
  mutate(., current_former_percentage = current_former_smoker/total*100) %>%
  mutate(., never_smoker_percentage = never_smoker/total*100) %>%
  mutate(., never_smoker_unknown_percentage = never_smoker_unknown/total*100) %>%
  mutate(., not_stated_percentage = not_stated/total*100) %>%
  mutate(., missing_percentage = missing/total*100) %>%
  select(lead_author, date_published, country, sample_size, median_age, iqr_lower, iqr_upper, mean_age, lower_range,
         upper_range, standard_deviation, female_sex_percent, current_percentage, former_percentage, current_former_percentage,
         never_smoker_percentage, never_smoker_unknown_percentage, not_stated_percentage, missing_percentage, study_id)

a <- data_study_general %>%
  select(study_id, study_setting)

table_1_word <- left_join(table_1_word, a, by = 'study_id') %>%
  select(1:4, 21, 5:20)
table_1_word$date_published <- as.Date.character(table_1_word$date_published)
write_rds(table_1_word, here::here('data_clean', 'table_1_word.rds'))

table_1_word <- table_1_word %>%
  mutate(median_mean = ifelse(is.na(median_age), mean_age, median_age))


a <- table_1_word %>%
  mutate(median_mean = ifelse(is.na(median_age), mean_age, median_age),
         mean_used = ifelse(is.na(mean_age), '','^'),
         iqr = ifelse(is.na(iqr_lower), NA, paste(iqr_lower, iqr_upper, sep = '-')),
         range_combined = paste(lower_range, upper_range, sep = '-'),
         range_combined = na_if(range_combined, 'NA-NA'),
         st_dev = paste((as.integer(median_mean-standard_deviation)), as.integer((median_mean+standard_deviation)), sep = '-'),
         st_dev = na_if(st_dev, 'NA-NA'))

a$iqr <- coalesce(a$iqr, a$range_combined, a$st_dev)
a$iqr_used <- coalesce(a$range_used, a$st_dev_used)

a <- a %>%
  select(study_id, lead_author, date_published, country, sample_size, study_setting, median_mean, mean_used, iqr,
         female_sex_percent, current_percentage, current_former_percentage, never_smoker_percentage,
         never_smoker_unknown_percentage, not_stated_percentage, missing_percentage) %>%
  mutate(median_mean = paste(median_mean, mean_used, sep = ''),
         median_mean = na_if(median_mean, NA)) %>%
  select(-mean_used) %>%
  mutate(median_mean = ifelse(median_mean == 'NA', 'NA', paste(paste(median_mean, iqr, sep = ' ('),')', sep = ''))) %>%
  select(-iqr) %>%
  rename('Study ID' = study_id,
         'Lead author' = lead_author,
         'Date published' = date_published,
         'Country' = country,
         'Sample size' = sample_size,
         'Study setting' = study_setting,
         'Median (IQR)' = median_mean,
         'Female %' = female_sex_percent,
         'Current smoker %' = current_percentage,
         'Current/former smokers %' = current_former_percentage,
         'Never smokers %' = never_smoker_percentage,
         'Never/unknown smokers %' = never_smoker_unknown_percentage,
         'Not stated %' = not_stated_percentage,
         'Missing %' = missing_percentage)
         
a$`Lead author` <- to_upper_camel_case(a$`Lead author`, sep_out = ", ")
a$`Lead author` <- a$`Lead author` %>%
  recode('Chow, Us, Cdc' = 'Chow (US CDC)',
         "Gold, Us, Cdc" = "Gold (US CDC)",
         "Miyara, Updated" = "Miyara",
         "De, La, Rica" = "de la Rica",
         "Opensafely, Collaborative" = "The Opensafely Collaborative",
         "Bello, Chavolla" = "Bello-Chavolla",
         "Carrillo, Vega" = "Carillo-Vega",
         "De, Lusignan" = "de Lusignan",
         "Mejia, Vilet" = "Mejia-Vilet",
         "Heili, Frades" = "Heili-Frades",
         "Vaquero, Roncero" = "Vaquero-Roncero",
         "Al, Hindawi" = "Al-Hindawi",
         "Del, Valle" = "del Valle",
         "Soto, Mota" = "Soto-Mota")

a$`Study setting` <-to_title_case(a$`Study setting`, sep_out = " ")
a$`Country` <-to_title_case(a$`Country`, sep_out = " ")
a$`Country` <- a$Country %>%
  recode('Usa' = 'USA',
         'Uk' = 'UK')

numeric_columns <- c('Median (IQR)', 'Female %', 'Current smoker %', 'Current/former smokers %',
                     'Never smokers %', 'Never/unknown smokers %', 'Not stated %', 'Missing %')

a <- flextable(a) %>%
    set_caption(caption = 'Characteristics of included studies') %>%
  colformat_num(col_keys = numeric_columns, digits = 2, na_str = '-', big.mark = ',') %>%
  colformat_num(col_keys = 'Sample size', digits = 0, na_str = '-', big.mark = ',') %>%
  set_table_properties(width = 1, layout = 'autofit')
a

Characteristics of included studies

Study ID	Lead author	Date published	Country	Sample size	Study setting	Median (IQR)	Female %	Current smoker %	Current/former smokers %	Never smokers %	Never/unknown smokers %	Not stated %	Missing %
1	Guan, Ni	2020-02-28	China	1,099	Hospital	47 (35-58)	41.90	12.47	-	84.35	-	-	1.27
2	Guan, Liang	2020-03-26	China	1,590	Hospital	49 (33-64)	42.70	-	6.98	93.02	-	-	-
3	Lian	2020-03-25	China	788	Hospital	NA	38.50	6.85	-	-	-	93.15	-
4	Jin	2020-03-24	China	651	Hospital	46 (32-60)	49.20	6.30	-	-	-	93.70	-
5	Chen	2020-03-26	China	548	Hospital	62 (44-70)	37.60	4.38	-	-	-	93.07	-
6	Zhou, Yu	2020-03-11	China	191	Hospital	56 (46-67)	38.00	5.76	-	-	-	94.24	-
7	Mo	2020-03-16	China	155	Hospital	54 (53-66)	44.50	3.87	-	-	-	96.13	-
8	Zhang, Dong	2020-02-19	China	140	Hospital	57^ (25-87)	46.30	1.43	-	-	-	93.57	-
9	Wan	2020-03-21	China	135	Hospital	47 (36-55)	46.70	6.67	-	-	-	93.33	-
10	Liu, Tao	2020-02-28	China	78	Hospital	38 (33-57)	50.00	-	6.41	-	-	93.59	-
11	Huang, Wang	2020-01-24	China	41	Hospital	49 (41-58)	27.00	7.32	-	-	-	92.68	-
12	Zhang, Cai	2020-03-20	China	645	Hospital	NA	49.10	6.36	-	-	-	93.64	-
13	Guo	2020-03-27	China	187	Hospital	59 (45-73)	51.30	9.63	-	-	-	90.37	-
14	Liu, Ming	2020-03-12	China	41	Hospital	39 (30-48)	58.50	9.76	-	-	-	90.24	-
15	Huang, Yang	2020-03-05	China	36	Hospital	69 (60-78)	30.60	-	11.11	-	-	88.89	-
16	Xu	2020-03-08	China	53	Hospital	NA	47.20	11.32	-	-	-	88.68	-
17	Li	2020-02-12	China	17	Hospital	45 (33-57)	47.10	17.65	-	-	-	82.35	-
18	Rentsch	2020-04-14	USA	3,528	Community and Hospital	66 (60-70)	4.60	27.18	-	36.92	-	-	5.30
19	Hu	2020-03-25	China	323	Hospital	61^ (23-91)	48.60	-	11.76	-	-	88.24	-
20	Wang, Pan	2020-03-24	China	125	Hospital	41 (26-66)	43.20	-	12.80	-	-	87.20	-
22	Petrilli	2020-04-11	USA	4,103	Community and Hospital	52 (36-65)	47.90	5.17	-	-	78.60	-	-
23	Chow (US CDC)	2020-03-31	USA	7,162	Community and Hospital	NA	-	1.34	-	-	-	-	96.36
25	Dong, Cao	2020-03-20	China	9	Hospital	44 (30-46)	66.70	11.11	-	-	-	88.89	-
26	Kim	2020-04-01	Korea	28	Hospital	43 (30-56)	46.40	17.86	-	-	-	78.57	3.57
27	Shi, Yu	2020-03-18	China	487	Hospital	46 (27-65)	46.80	-	8.21	-	-	89.12	2.67
28	Yang, Yu	2020-02-24	China	52	Hospital	60 (47-73)	37.00	3.85	-	-	-	96.15	-
29	Argenziano	2020-04-22	USA	1,000	Hospital	63 (50-75)	40.40	4.90	-	77.20	-	-	-
30	Solis	2020-04-25	Mexico	650	Hospital	46 (NA)	42.10	9.38	-	-	-	90.62	-
31	Richardson	2020-04-22	USA	5,700	Hospital	63 (52-75)	39.70	-	9.79	52.79	-	-	37.42
32	Fontanet	2020-04-23	France	661	Community and Hospital	37 (16-47)	62.00	10.44	-	-	89.56	-	-
33	Zheng, Gao	2020-04-19	China	66	Hospital	47^ (NA)	25.80	12.12	-	-	-	87.88	-
34	Liao, Feng	2020-04-24	China	1,848	Hospital	55 (48-61)	54.70	-	0.43	-	-	7.58	91.99
35	Rodriguez	2020-04-24	Spain	7	Hospital	68 (34-75)	28.60	-	42.86	57.14	-	-	-
36	Magagnoli	2020-04-16	USA	368	Hospital	69 (59-75)	0.00	-	14.13	-	-	85.87	-
37	Shi, Ren	2020-04-23	China	134	Hospital	46 (34-58)	51.50	-	10.45	-	-	89.55	-
38	Hadjadj	2020-04-23	France	50	Hospital	55 (50-63)	22.00	2.00	-	80.00	-	-	-
39	Niedzwiedz	2020-04-30	UK	1,474	Community and Hospital	NA	-	9.98	-	55.04	-	-	0.59
41	Gold (US CDC)	2020-04-20	USA	305	Hospital	NA	50.50	5.25	-	-	-	-	94.75
44	Yu, Cai	2020-04-27	China	95	Hospital	NA	44.21	8.42	-	-	-	91.58	-
45	Zheng, Xiong	2020-04-30	China	73	Hospital	43^ (NA)	45.20	-	10.96	89.04	-	-	-
46	Miyara	2020-05-09	France	479	Community and Hospital	NA	44.70	6.68	-	59.71	-	-	1.88
47	de la Rica	2020-05-11	Spain	48	Hospital	66^ (33-88)	33.00	-	20.83	-	-	77.08	2.08
48	Yin, Yang	2020-05-10	China	106	Hospital	73 (61-85)	39.60	-	16.98	-	-	83.02	-
49	Gaibazzi	2020-05-10	Italy	441	Hospital	71 (62-80)	38.00	4.76	-	85.26	-	-	-
50	Shi, Zuo	2020-05-10	USA	96	Hospital	63^ (44-82)	41.00	-	30.21	-	-	69.79	-
51	Cho	2020-05-11	UK	1,331	Community and Hospital	NA	49.20	19.01	-	54.02	-	-	-
52	Allenbach	2020-05-08	France	152	Hospital	77 (60-83)	31.10	-	6.58	-	-	93.42	-
53	Robilotti	2020-05-08	USA	423	Hospital	NA	50.00	2.13	-	58.63	-	-	1.65
54	The Opensafely Collaborative	2020-05-07	UK	17,425,445	Community and Hospital	NA	50.10	17.00	-	45.91	-	-	4.16
55	Borobia	2020-05-06	Spain	2,226	Hospital	61 (46-78)	52.00	7.05	-	-	-	92.95	-
56	Giacomelli	2020-05-06	Italy	233	Hospital	61 (50-72)	31.90	-	30.04	69.96	-	-	-
57	Shah	2020-05-06	USA	316	Hospital	63 (43-72)	48.10	16.46	-	42.09	-	-	23.73
58	Bello-Chavolla	2020-05-06	Mexico	62,489	Community and Hospital	NA	49.40	-	9.94	-	-	90.06	-
59	Kolin	2020-05-05	UK	1,474	Community and Hospital	58 (49-67)	46.60	14.45	-	44.57	-	-	0.81
60	Lubetzky	2020-05-08	USA	54	Hospital	57 (29-83)	62.00	-	22.22	-	-	77.78	-
61	Goyal	2020-04-17	USA	393	Hospital	62.2 (49-74)	39.30	5.09	-	-	-	94.91	-
62	Feng	2020-04-10	China	476	Hospital	53 (40-64)	43.10	9.24	-	-	-	86.13	4.62
63	Yao	2020-04-24	China	108	Hospital	52 (37-58)	60.20	3.70	-	-	-	-	96.30
64	Sami	2020-05-15	Iran	490	Hospital	56.6 (41-71)	39.00	14.08	-	-	85.92	-	-
65	Almazeedi	2020-05-15	Kuwait	1,096	Hospital	41 (25-57)	19.00	4.01	-	-	95.99	-	-
66	Carillo-Vega	2020-05-14	Mexico	10,544	Community and Hospital	46.5^ (30-62)	42.30	8.88	-	-	-	91.12	-
67	Yanover	2020-05-13	Israel	4,353	Community and Hospital	35 (22-54)	44.50	11.81	-	85.23	-	-	-
68	Hamer	2020-05-13	UK	387,109	Hospital	56.2 (48-64)	55.10	9.67	-	55.50	-	-	-
69	Regina	2020-05-14	Switzerland	200	Hospital	70 (55-81)	40.00	4.50	-	-	-	95.50	-
71	de Lusignan	2020-05-15	UK	3,802	Community	58 (34-73)	57.60	10.86	-	29.59	-	-	13.44
72	Targher	2020-05-13	China	339	Hospital	48.4^ (NA)	52.80	8.26	-	-	-	91.74	-
73	Valenti	2020-05-18	Italy	789	Community	40.7^ (NA)	34.98	25.86	-	-	-	74.14	-
74	Feuth	2020-05-18	Finland	28	Hospital	56 (47-72)	46.00	10.71	-	60.71	-	-	-
75	Heng	2020-05-18	China	51	Hospital	70 (58-79)	27.50	13.73	-	-	-	86.27	-
76	Parrotta	2020-05-18	USA	76	Community and Hospital	44.9 (13-71)	61.80	2.63	-	68.42	-	2.63	-
77	Shekar	2020-05-18	USA	50	Hospital	55.5 (20-85)	54.00	48.00	-	-	-	52.00	-
78	Mejia-Vilet	2020-05-16	Mexico	329	Hospital	49 (41-60)	36.00	-	6.99	-	-	93.01	-
79	Chen, Jiang	2020-05-16	China	135	Hospital	NA	42.20	-	9.63	-	-	90.37	-
80	Jiong	2020-05-16	China	1,008	Hospital	55 (44-65)	43.60	5.65	-	-	-	48.02	46.33
81	Rimland	2020-05-19	USA	11	Hospital	59 (48-65)	18.20	9.09	-	-	-	81.82	-
82	Palaiodimos	2020-05-15	USA	200	Hospital	64 (50-73.5)	51.00	-	32.50	67.50	-	-	-
83	Ip	2020-05-25	USA	2,512	Hospital	64 (52-76)	37.62	3.07	-	64.49	-	-	14.61
84	Heili-Frades	2020-05-25	Spain	4,712	Hospital	62 (47-77)	50.50	4.94	-	-	66.49	-	11.16
85	Vaquero-Roncero	2020-05-24	Spain	146	Hospital	66^ (59-72)	32.20	-	6.85	-	-	93.15	-
86	Kim, Garg	2020-05-22	USA	2,491	Hospital	62 (50-75)	46.80	6.02	-	-	68.13	-	0.08
87	Wu	2020-05-21	Italy	174	Hospital	61.2^ (50-71)	30.46	-	33.33	-	-	66.67	-
88	Hopkinson	2020-05-21	UK	2,401,982	Community	43.6^ (32-55)	63.30	-	10.96	-	-	89.04	-
89	Shi, Zhao	2020-05-20	China	101	Hospital	71 (59-80)	40.60	-	4.95	-	-	95.05	-
90	Kimmig	2020-05-20	USA	60	Hospital	64 (50-68)	41.67	-	36.67	-	-	63.33	-
91	Al-Hindawi	2020-05-20	UK	31	Hospital	61 (NA)	12.90	3.23	-	25.81	-	-	-
92	Basse	2020-05-19	France	141	Hospital	62 (52-72)	72.00	17.73	-	-	-	82.27	-
93	Freites	2020-05-19	Spain	123	Hospital	59.88^ (44-74)	69.92	3.25	-	-	-	96.75	-
94	Alshami	2020-05-19	Saudi Arabia	128	Quarantine Centre	39.6^ (24-55)	53.90	15.62	-	-	-	82.03	-
95	Russell	2020-05-19	UK	106	Hospital	67^ (52-81)	45.30	6.60	-	39.62	-	29.25	-
96	Berumen	2020-05-26	Mexico	102,875	Hospital	NA	49.08	-	9.64	-	90.36	-	-
97	Gianfrancesco	2020-05-29	Multiple	600	Community and Hospital	56 (45-67)	71.00	-	21.50	64.83	-	-	13.67
98	Li, Long	2020-05-28	China	145	Not Stated	49^ (13-80)	61.00	-	5.52	-	-	94.48	-
99	Batty	2020-06-01	UK	908	Hospital	57.27^ (48-66)	44.27	11.23	-	-	-	88.77	-
100	Israel	2020-06-01	Israel	24,087	Community and Hospital	43.4^ (24-62)	48.70	17.08	-	69.88	-	-	-
101	del Valle	2020-05-30	USA	1,484	Hospital	62 (52-72)	40.60	5.53	-	-	-	45.82	25.34
102	Zuo, Zuo	2020-05-29	USA	44	Hospital	57^ (45-69)	18.18	-	27.27	-	-	72.73	-
103	Chaudhry	2020-05-29	USA	40	Community and Hospital	52 (45.5-61)	60.00	-	15.00	-	-	82.50	2.50
104	Louis	2020-06-28	USA	22	Hospital	66.5^ (55-77)	36.40	-	45.45	-	-	54.55	-
105	Soto-Mota	2020-05-27	Mexico	400	Hospital	NA	30.00	-	12.00	-	-	88.00	-
106	Patel	2020-05-26	USA	104	Hospital	60.66^ (47-74)	47.00	41.35	-	-	49.04	-	9.62
107	Garibaldi	2020-05-26	USA	832	Hospital	63 (49-75)	47.00	5.53	-	-	-	71.88	-
108	Docherty	2020-05-22	Multiple	20,133	Hospital	72.9 (58-82)	40.00	4.23	-	44.54	-	-	29.55
109	Boulware	2020-03-06	Multiple	821	Community	40 (33-50)	51.60	3.29	-	-	-	96.71	-

save_as_docx(a, path = here('data_clean', 'Table_1.docx'))

0.5.2 Table 2: Smoking and COVID-19 diagnostic testing

quality_rating <- table_6 %>%
  select(author, 16:17) %>%
  rename('lead_author' = author, 'overall_rating' = 2, 'missingness' = 3) %>%
  clean_data()

#Table 2
table_2_word <-  table_2 %>%
  mutate(., sample = contributing_sample) %>%
  mutate(., negative_test_percentage =formatC(negative_test/sample*100, digits = 2, format = "f"),
         negative_current_percentage = formatC(negative_current_smoker/negative_test*100, digits = 2, format = "f"),
         negative_former_smoker_percentage = formatC(negative_former_smoker/negative_test*100, digits = 2, format = "f"),
         negative_current_former_smoker_percentage = formatC(negative_current_former_smoker/negative_test*100, digits = 2, format = "f"),
         negative_never_smoker_percentage = formatC(negative_never_smoker/negative_test*100, digits = 2, format = "f"),
         negative_not_stated_percentage = formatC(negative_not_stated/negative_test*100, digits = 2, format = "f"),
         positive_test_percentage = formatC(positive_test/sample*100, digits = 2, format = "f"),
         positive_current_smoker_percentage = formatC(positive_current_smoker/positive_test*100, digits = 2, format = "f"),
         positive_former_smoker_percentage = formatC(positive_former_smoker/positive_test*100, digits = 2, format = "f"),
         positive_current_former_smoker_percentage = formatC(positive_current_former_smoker/positive_test*100, digits = 2, format = "f"),
         positive_never_smoker_percentage = formatC(positive_never_smoker/positive_test*100, digits = 2, format = "f"),
         positive_not_stated_percentage = formatC(positive_not_stated/positive_test*100, digits = 2, format = "f")) %>%
  select(-data_on_testing, -missing, -date_published, -sample) %>%
  write_rds(here::here('data_clean', 'table_2_word.rds'))

quality_table_2 <- table_2_word %>%
  left_join(., quality_rating, by = 'lead_author') %>%
  select(lead_author, overall_rating, missingness) %>%
  mutate(overall_quality = ifelse(overall_rating != 'poor' & missingness < 20, 'fair', 'poor')) %>%
  select(lead_author, overall_quality)

a <- table_2_word %>%
  filter(contributing_sample >= 1) %>%
  mutate(Author = lead_author,
         Population_tested = contributing_sample,
         SARS_CoV_2_negative = paste(
           paste(negative_test, 
                 (negative_test_percentage), sep = " ("), "%)", sep = ""),
         N_current_smoker = paste(paste(negative_current_smoker, (negative_current_percentage), sep = " (")
                                  , "%)", sep = ""),
         N_current_smoker = na_if(N_current_smoker, "NA ( NA%)"),
         N_former_smoker = paste(paste(negative_former_smoker, (negative_former_smoker_percentage), sep = " (")
                                 , "%)", sep = ""),
         N_former_smoker = na_if(N_former_smoker, "NA ( NA%)"),
         N_current_former_smoker = paste(paste(negative_current_former_smoker,
                                               (negative_current_former_smoker_percentage), sep = " (")
                                         , "%)", sep = ""),
         N_current_former_smoker = na_if(N_current_former_smoker, "NA ( NA%)"),
         N_never_smoker = paste(paste(negative_never_smoker,
                                      (negative_never_smoker_percentage), sep = " (")
                                , "%)", sep = ""),
         N_never_smoker = na_if(N_never_smoker, "NA ( NA%)"),
         N_not_stated = paste(paste(negative_not_stated,
                                    (negative_not_stated_percentage), sep = " (")
                              , "%)", sep = ""),
         N_not_stated = na_if(N_not_stated, "NA ( NA%)")) %>%
  mutate(SARS_CoV_2_positive = paste(
    paste(positive_test, 
          (positive_test_percentage), sep = " ("), "%)", sep = ""),
    P_current_smoker = paste(paste(positive_current_smoker, (positive_current_smoker_percentage), sep = " (")
                             , "%)", sep = ""),
    P_current_smoker = na_if(P_current_smoker, "NA ( NA%)"),
    P_former_smoker = paste(paste(positive_former_smoker, (positive_former_smoker_percentage), sep = " (")
                            , "%)", sep = ""),
    P_former_smoker = na_if(P_former_smoker, "NA ( NA%)"),
    P_current_former_smoker = paste(paste(positive_current_former_smoker,
                                          (positive_current_former_smoker_percentage), sep = " (")
                                    , "%)", sep = ""),
    P_current_former_smoker = na_if(P_current_former_smoker, "NA ( NA%)"),
    P_never_smoker = paste(paste(positive_never_smoker,
                                 (positive_never_smoker_percentage), sep = " (")
                           , "%)", sep = ""),
    P_never_smoker = na_if(P_never_smoker, "NA ( NA%)"),
    P_not_stated = paste(paste(positive_not_stated,
                               (positive_not_stated_percentage), sep = " (")
                         , "%)", sep = ""),
    P_not_stated = na_if(P_not_stated, "NA ( NA%)")) %>%
  select(Author, Population_tested, SARS_CoV_2_negative, N_current_smoker, N_former_smoker, N_current_former_smoker,
         N_never_smoker, N_not_stated, SARS_CoV_2_positive, P_current_smoker, P_former_smoker, P_current_former_smoker,
         P_never_smoker, P_not_stated)

         
a$Author <- to_upper_camel_case(a$Author, sep_out = ", ")
a$Author  <- a$Author  %>%
  recode("Bello, Chavolla" = "Bello-Chavolla",
         "De, Lusignan" = "de Lusignan",
         "Del, Valle" = "del Valle")

numeric_columns <- c('Population_tested', 'SARS_CoV_2_negative', 'N_current_smoker', 'N_former_smoker',
                     'N_current_former_smoker', 'N_never_smoker', 'N_not_stated', 'SARS_CoV_2_positive', 'P_current_smoker',
                     'P_former_smoker', 'P_current_former_smoker', 'P_never_smoker', 'P_not_stated')
a <- flextable(a) %>%
    set_caption(caption = 'SARS-CoV-2 infection by smoking status') %>%
  colformat_num(col_keys = numeric_columns, digits = 0, na_str = '-', big.mark = ',')

a <- set_header_labels(a,
                       Population_tested = 'Total population tested',
                       SARS_CoV_2_negative = "N (%)",
                       N_current_smoker = "Current smoker (%)",
                       N_former_smoker = "Former smoker (%)",
                       N_current_former_smoker = "Current/former smoker (%)",
                       N_never_smoker = "Never smoker (%)",
                       N_not_stated = "Not stated (%)",
                       SARS_CoV_2_positive = "N (%)",
                       P_current_smoker = "Current smoker (%)",
                       P_former_smoker = "Former smoker (%)",
                       P_current_former_smoker = "Current/former smoker (%)",
                       P_never_smoker = "Never smoker (%)",
                       P_not_stated = "Not stated (%)") %>%
  add_header_row(top = TRUE, values = c("","SARS-CoV-2 negative", "SARS-CoV-2 positive" ), colwidths = c(2, 6, 6)) %>%
  theme_booktabs() %>%
  fix_border_issues() %>%
  set_table_properties(width = 1, layout = 'autofit')  %>%
  colformat_char(na_str = '-')
a

SARS-CoV-2 infection by smoking status

		SARS-CoV-2 negative						SARS-CoV-2 positive
Author	Total population tested	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Not stated (%)	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Not stated (%)
Rentsch	3528	2974 (84.30%)	1444 (48.55%)	704 (23.67%)	-	826 (27.77%)	-	554 (15.70%)	159 (28.70%)	179 (32.31%)	-	216 (38.99%)	-
Fontanet	661	490 (74.13%)	64 (13.06%)	-	-	426 (86.94%)	-	171 (25.87%)	5 (2.92%)	-	-	166 (97.08%)	-
Cho	1331	793 (59.58%)	142 (17.91%)	214 (26.99%)	-	437 (55.11%)	-	538 (40.42%)	111 (20.63%)	145 (26.95%)	-	282 (52.42%)	-
Shah	243	212 (87.24%)	52 (24.53%)	47 (22.17%)	-	113 (53.30%)	-	29 (11.93%)	0 (0.00%)	9 (31.03%)	-	20 (68.97%)	-
Bello-Chavolla	62489	46960 (75.15%)	-	-	4835 (10.30%)	-	42125 (89.70%)	15529 (24.85%)	-	-	1374 (8.85%)	-	14155 (91.15%)
Kolin	1474	805 (54.61%)	141 (17.52%)	307 (38.14%)	-	354 (43.98%)	3 (0.37%)	669 (45.39%)	72 (10.76%)	285 (42.60%)	-	303 (45.29%)	9 (1.35%)
de Lusignan	3291	2740 (83.26%)	366 (13.36%)	1450 (52.92%)	-	924 (33.72%)	-	551 (16.74%)	47 (8.53%)	303 (54.99%)	-	201 (36.48%)	-
Valenti	789	689 (87.33%)	197 (28.59%)	-	-	-	492 (71.41%)	40 (5.07%)	7 (17.50%)	-	-	-	33 (82.50%)
Parrotta	76	39 (51.32%)	1 (2.56%)	10 (25.64%)	-	27 (69.23%)	1 (2.56%)	37 (48.68%)	1 (2.70%)	10 (27.03%)	-	25 (67.57%)	1 (2.70%)
Hopkinson	26918	16759 (62.26%)	-	-	1552 (9.26%)	-	15207 (90.74%)	7123 (26.46%)	-	-	527 (7.40%)	-	6596 (92.60%)
Berumen	102875	71353 (69.36%)	-	-	7173 (10.05%)	64180 (89.95%)	-	31522 (30.64%)	-	-	2748 (8.72%)	28774 (91.28%)	-
Israel	24087	20076 (83.35%)	3711 (18.48%)	2670 (13.30%)	-	13695 (68.22%)	-	4011 (16.65%)	403 (10.05%)	471 (11.74%)	-	3137 (78.21%)	-
del Valle	1108	143 (12.91%)	27 (18.88%)	53 (37.06%)	-	-	63 (44.06%)	965 (87.09%)	55 (5.70%)	293 (30.36%)	-	-	617 (63.94%)

nrow(table_2_word)

## [1] 14

save_as_docx(a, path = here('data_clean', 'Table_2.docx'))

0.5.2.1 Smoking and testing meta-analysis

source(here::here('scripts', 'rr_function.R'))

table_2 <- table_2_word

included_studies <- quality_table_2 %>%
  filter(overall_quality == 'fair') %>%
  select(lead_author)

meta <- tibble('author' = table_2$lead_author,
               'negative_smoker' = table_2$negative_current_smoker,
               'negative_never_smoker' = table_2$negative_never_smoker,
               'positive_smoker' = table_2$positive_current_smoker,
               'positive_never_smoker' = table_2$positive_never_smoker,
               'negative_former_smoker' = table_2$negative_former_smoker,
               'positive_former_smoker' = table_2$positive_former_smoker) %>%
        filter(author %in% included_studies$lead_author)

meta$author <- recode(meta$author, 'rentsch' = 'Rentsch',
                      'cho' = 'Cho',
                      'kolin' = 'Kolin',
                      'de_lusignan' = 'de Lusignan',
                      'parrotta' = 'Parrotta',
                      'israel' = 'Israel')

a <- RR_testing('Rentsch', 'current')
b <- RR_testing('Rentsch', 'former')
c <- RR_testing('Cho', 'current')
d <- RR_testing('Cho', 'former')
e <- RR_testing('Kolin', 'current')
f <- RR_testing('Kolin', 'former')
g <- RR_testing('de Lusignan', 'current')
h <- RR_testing('de Lusignan', 'former')
i <- RR_testing('Parrotta', 'current')
j <- RR_testing('Parrotta', 'former')
k <- RR_testing('Israel', 'current')
l <- RR_testing('Israel', 'former')

data <- list(a,b,c,d,e,f,g,h,i,j,k,l)
k <- do.call(rbind.data.frame, data)

#SEs for Niedzwiedz et al. 2020

#current vs. never

niedz_log_RR_1<-log(1.15)
niedz_log_SE_1<-(log(1.54)-log(0.86))/3.92

k <- k %>%
  add_row(., study = 'Niedzwiedz', smoking_status = 'current', log_RR = niedz_log_RR_1, log_SE = niedz_log_SE_1)
#former vs. never

niedz_log_RR_2<-log(1.42)
niedz_log_SE_2<-(log(1.69)-log(1.19))/3.92

k <- k %>%
  add_row(., study = 'Niedzwiedz', smoking_status = 'former', log_RR = niedz_log_RR_2, log_SE = niedz_log_SE_2)

numbers_in_analysis <- table_2_word %>%
        left_join(., quality_rating, by = 'lead_author') %>%
        filter(., overall_rating != 'poor') %>%
        add_row(lead_author = 'niedzwiedz', contributing_sample = 1474) %>%
        replace_na(list(negative_not_stated = 0, positive_not_stated = 0)) %>%
        mutate(contributing_sample = (contributing_sample - (negative_not_stated+positive_not_stated))) %>%
                       select(lead_author, contributing_sample) %>%
  na.omit()

The studies included in meta-analysis are c(“rentsch”, “niedzwiedz”, “cho”, “kolin”, “de_lusignan”, “parrotta”, “israel”).

• Current versus never smokers

#current vs. never smokers
current_never_meta <- k %>%
  filter(smoking_status == 'current')
  
a<-metagen(current_never_meta$log_RR,
           current_never_meta$log_SE,
           studlab = c("Rentsch, n = 3528",
                       "Cho, n = 1331", 
                       "Kolin, = 1462", 
                       "de Lusignan, n = 3291", 
                       "Parrotta, n = 74",
                       "Israel, n = 24087 ", 
                       "Niedzwiedz, n = 1474"),
           sm="RR",
           comb.fixed = F, comb.random = T)

png(here::here('reports', 'figure','fig_3.png'), width=1024, height=546, res=120)
forest(a,
       xlim = c(0.3, 2),
       sortvar = a$TE,
       rightlabs = c('RR', '95% CI', 'Weight'),
       leftlabs = c('Author', 'logRR', 'SE'),
       print.tau2 = F,
       col.diamond = 'blue',
       col.diamond.lines = 'black',
       col.square = 'black',
       col.square.lines = 'black',
       digits.sd = 2,
       colgap.forest.left = unit(15,"mm"))
null <- dev.off()

a

##                           RR           95%-CI %W(random)
## Rentsch, n = 3528     0.4785 [0.3959; 0.5783]       16.4
## Cho, n = 1331         1.1186 [0.9471; 1.3212]       16.6
## Kolin, = 1462         0.7330 [0.5969; 0.9000]       16.2
## de Lusignan, n = 3291 0.6369 [0.4733; 0.8571]       15.0
## Parrotta, n = 74      1.0400 [0.2528; 4.2784]        3.6
## Israel, n = 24087     0.5256 [0.4766; 0.5797]       17.2
## Niedzwiedz, n = 1474  1.1500 [0.8594; 1.5389]       15.1
## 
## Number of studies combined: k = 7
## 
##                          RR           95%-CI     z p-value
## Random effects model 0.7350 [0.5430; 0.9948] -1.99  0.0462
## 
## Quantifying heterogeneity:
##  tau^2 = 0.1363 [0.0342; 0.7515]; tau = 0.3692 [0.1849; 0.8669];
##  I^2 = 93.0% [88.0%; 95.9%]; H = 3.77 [2.89; 4.93]
## 
## Test of heterogeneity:
##      Q d.f.  p-value
##  85.44    6 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
## - Jackson method for confidence interval of tau^2 and tau

Current versus never smokers

• Former versus never smokers

#former vs. never smokers

current_never_meta <- k %>%
  filter(smoking_status == 'former')

a<-metagen(current_never_meta$log_RR,
           current_never_meta$log_SE,
           studlab = c("Rentsch, n = 3528",
                       "Cho, n = 1331", 
                       "Kolin, = 1462", 
                       "de Lusignan, n = 3291", 
                       "Parrotta, n = 74",
                       "Israel, n = 24087 ", 
                       "Niedzwiedz, n = 1474"),
           sm="RR", comb.fixed = F, comb.random = T)

png(here::here('reports', 'figure','fig_4.png'), width=1024, height=546, res=120)
forest(a,
       xlim = c(0.3, 2),
       sortvar = a$TE,
       rightlabs = c('RR', '95% CI', 'Weight'),
       leftlabs = c('Author', 'logRR', 'SE'),
       print.tau2 = F,
       col.diamond = 'blue',
       col.diamond.lines = 'black',
       col.square = 'black',
       col.square.lines = 'black',
       digits.sd = 2,
       colgap.forest.left = unit(15,"mm"))
null <- dev.off()

a

##                           RR           95%-CI %W(random)
## Rentsch, n = 3528     0.9779 [0.8196; 1.1669]       14.7
## Cho, n = 1331         1.0298 [0.8818; 1.2026]       15.5
## Kolin, = 1462         1.0439 [0.9281; 1.1741]       16.7
## de Lusignan, n = 3291 0.9674 [0.8229; 1.1373]       15.2
## Parrotta, n = 74      1.0400 [0.6174; 1.7517]        5.6
## Israel, n = 24087     0.8046 [0.7360; 0.8795]       17.5
## Niedzwiedz, n = 1474  1.4200 [1.1916; 1.6922]       14.8
## 
## Number of studies combined: k = 7
## 
##                          RR           95%-CI    z p-value
## Random effects model 1.0195 [0.8803; 1.1807] 0.26  0.7968
## 
## Quantifying heterogeneity:
##  tau^2 = 0.0300 [0.0048; 0.1616]; tau = 0.1732 [0.0693; 0.4020];
##  I^2 = 83.8% [68.3%; 91.8%]; H = 2.49 [1.78; 3.48]
## 
## Test of heterogeneity:
##      Q d.f.  p-value
##  37.12    6 < 0.0001
## 
## Details on meta-analytical method:
## - Inverse variance method
## - DerSimonian-Laird estimator for tau^2
## - Jackson method for confidence interval of tau^2 and tau

Former versus never smokers

0.5.3 Table 3: Smoking and COVID-19 hospitalisation

#Table 3
table_3_word <- table_3 %>%
  mutate(., sample = sample_with_outcome) %>%
  mutate(., community_percentage = formatC(number_community/sample*100, digits = 2, format = "f")) %>%
  mutate(., community_current_smoker_percent = formatC(community_current_smoker/number_community*100, digits = 2, format = "f")) %>%
  mutate(., community_former_smoker_percent = formatC(community_former_smoker/number_community*100, digits = 2, format = "f")) %>%
  mutate(., community_current_former_smoker_percent = formatC(community_current_former_smoker/number_community*100, digits = 2, format = "f")) %>%
  mutate(., community_never_smoker_percent = formatC(community_never_smoker/number_community*100, digits = 2, format = "f")) %>%
  mutate(., community_never_unknown_smoker_percent = formatC(community_never_unknown_smoker/number_community*100, digits = 2, format = "f")) %>%
  mutate(., community_not_stated_percent = formatC(community_not_stated/number_community*100, digits = 2, format = "f")) %>%
  mutate(., number_hospitalised_percent = formatC(number_hospitalised/sample*100, digits = 2, format = "f")) %>%
  mutate(., hospitalised_current_smoker_percent = formatC(hospitalised_current_smoker/number_hospitalised*100, digits = 2, format = "f")) %>%
  mutate(., hospitalised_former_smoker_percent = formatC(hospitalised_former_smoker/number_hospitalised*100, digits = 2, format = "f")) %>%
  mutate(., hospitalised_current_former_smoker_percent = formatC(hospitalised_current_former_smoker/number_hospitalised*100, digits = 2, format = "f")) %>%
  mutate(., hospitalised_never_smoker_percent = formatC(hospitalised_never_smoker/number_hospitalised*100, digits = 2, format = "f")) %>%
  mutate(., hospitalised_never_unknown_smoker_percent = formatC(hospitalised_never_unknown_smoker/number_hospitalised*100, digits = 2, format = "f")) %>%
  mutate(., hospitalised_not_stated_percent = formatC(hospitalised_not_stated/number_hospitalised*100, digits = 2, format = "f")) %>%
  select(lead_author, sample_with_outcome, number_community, community_percentage,
         community_current_smoker, community_current_smoker_percent, community_former_smoker, 
         community_former_smoker_percent, community_current_former_smoker,
         community_current_former_smoker_percent, community_never_smoker, community_never_smoker_percent,
         community_never_unknown_smoker, community_never_unknown_smoker_percent, 
         community_not_stated, community_not_stated_percent, number_hospitalised, number_hospitalised_percent,
         hospitalised_current_smoker, hospitalised_current_smoker_percent, hospitalised_former_smoker,
         hospitalised_former_smoker_percent, hospitalised_current_former_smoker, hospitalised_current_former_smoker_percent,
         hospitalised_never_smoker, hospitalised_never_smoker_percent, hospitalised_never_unknown_smoker,
         hospitalised_never_unknown_smoker_percent, hospitalised_not_stated, hospitalised_not_stated_percent) %>%
  write_rds(here::here('data_clean', 'table_3_word.rds'))

quality_table_3 <- table_3_word %>%
  left_join(., quality_rating, by = 'lead_author') %>%
  select(lead_author, overall_rating, missingness) %>%
  mutate(overall_quality = ifelse(overall_rating != 'poor' & missingness < 20, 'fair', 'poor')) %>%
  select(lead_author, overall_quality)

a <- table_3_word %>%
  filter(number_community >= 1) %>%
  mutate(Author = lead_author,
         Population = sample_with_outcome,
         Community = paste(
           paste(number_community, 
                 as.integer(community_percentage), sep = " ("), "%)", sep = ""),
         C_current_smoker = paste(paste(community_current_smoker, (community_current_smoker_percent), sep = " (")
                                  , "%)", sep = ""),
         C_current_smoker = na_if(C_current_smoker, "NA ( NA%)"),
         C_former_smoker = paste(paste(community_former_smoker, (community_former_smoker_percent), sep = " (")
                                 , "%)", sep = ""),
         C_former_smoker = na_if(C_former_smoker, "NA ( NA%)"),
         C_current_former_smoker = paste(paste(community_current_former_smoker,
                                               (community_current_former_smoker_percent), sep = " (")
                                         , "%)", sep = ""),
         C_current_former_smoker = na_if(C_current_former_smoker, "NA ( NA%)"),
         C_never_smoker = paste(paste(community_never_smoker,
                                      (community_never_smoker_percent), sep = " (")
                                , "%)", sep = ""),
         C_never_smoker = na_if(C_never_smoker, "NA ( NA%)"),
         C_never_unknown_smoker = paste(paste(community_never_unknown_smoker,
                                              (community_never_unknown_smoker_percent), sep = " (")
                                        , "%)", sep = ""),
         C_never_unknown_smoker = na_if(C_never_unknown_smoker, "NA ( NA%)"),
         C_not_stated = paste(paste(community_not_stated,
                                    (community_not_stated_percent), sep = " (")
                              , "%)", sep = ""),
         C_not_stated = na_if(C_not_stated, "NA ( NA%)")) %>%
  mutate(Hospitalised = paste(
    paste(number_hospitalised, 
          as.integer(number_hospitalised_percent), sep = " ("), "%)", sep = ""),
    H_current_smoker = paste(paste(hospitalised_current_smoker, (hospitalised_current_smoker_percent), sep = " (")
                             , "%)", sep = ""),
    H_current_smoker = na_if(H_current_smoker, "NA ( NA%)"),
    H_former_smoker = paste(paste(hospitalised_former_smoker, (hospitalised_former_smoker_percent), sep = " (")
                            , "%)", sep = ""),
    H_former_smoker = na_if(H_former_smoker, "NA ( NA%)"),
    H_current_former_smoker = paste(paste(hospitalised_current_former_smoker,
                                          (hospitalised_current_former_smoker_percent), sep = " (")
                                    , "%)", sep = ""),
    H_current_former_smoker = na_if(H_current_former_smoker, "NA ( NA%)"),
    H_never_smoker = paste(paste(hospitalised_never_smoker,
                                 (hospitalised_never_smoker_percent), sep = " (")
                           , "%)", sep = ""),
    H_never_smoker = na_if(H_never_smoker, "NA ( NA%)"),
    H_never_unknown_smoker = paste(paste(hospitalised_never_unknown_smoker,
                                         (hospitalised_never_unknown_smoker_percent), sep = " (")
                                   , "%)", sep = ""),
    H_never_unknown_smoker = na_if(H_never_unknown_smoker, "NA ( NA%)"),
    H_not_stated = paste(paste(hospitalised_not_stated,
                               (hospitalised_not_stated_percent), sep = " (")
                         , "%)", sep = ""),
    H_not_stated = na_if(H_not_stated, "NA ( NA%)")) %>%
  select(Author, Population, Community, C_current_smoker, C_former_smoker, C_current_former_smoker,
         C_never_smoker,C_never_unknown_smoker,C_not_stated, Hospitalised, H_current_smoker, H_former_smoker,
         H_current_former_smoker, H_never_smoker, H_never_unknown_smoker, H_not_stated)

         
a$Author <- to_upper_camel_case(a$Author, sep_out = ", ")
a$Author  <- a$Author  %>%
  recode("Chow, Us, Cdc" = "Chow (US CDC)",
         "Miyara, Updated" = "Miyara",
         "Carrillo, Vega" = "Carrillo-Vega",
         "Heili, Frades" = "Heili-Frades")

numeric_columns <- 'Population'

a <- flextable(a) %>%
    set_caption(caption = 'COVID-19 hospitalisation by smoking status') %>%
  colformat_num(col_keys = numeric_columns, digits = 0, na_str = '-', big.mark = ',') %>%
  colformat_char(na_str = '-')

a <- set_header_labels(a,
                       Population = 'Population with outcome',
                       Community = "N (%)",
                       C_current_smoker = "Current smoker (%)",
                       C_former_smoker = "Former smoker (%)",
                       C_current_former_smoker = "Current/former smoker (%)",
                       C_never_smoker = "Never smoker (%)",
                       C_never_unknown_smoker = "Never/unknown smoker (%)",
                       C_not_stated = "Not stated (%)",
                       Hospitalised = "N (%)",
                       H_current_smoker = "Current smoker (%)",
                       H_former_smoker = "Former smoker (%)",
                       H_current_former_smoker = "Current/former smoker (%)",
                       H_never_smoker = "Never smoker (%)",
                       H_never_unknown_smoker = "Never/unknown smoker (%)",
                       H_not_stated = "Not stated (%)") %>%
  add_header_row(top = TRUE, values = c("","Community", "Hospitalised" ), colwidths = c(2, 7, 7)) %>%
  theme_booktabs() %>%
  fix_border_issues() %>%
  set_table_properties(width = 1, layout = 'autofit')
a

COVID-19 hospitalisation by smoking status

		Community							Hospitalised
Author	Population with outcome	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Never/unknown smoker (%)	Not stated (%)	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Never/unknown smoker (%)	Not stated (%)
Rentsch	554	269 (48%)	69 (25.65%)	90 (33.46%)	-	110 (40.89%)	-	-	285 (51%)	90 (31.58%)	89 (31.23%)	-	106 (37.19%)	-	-
Petrilli	4103	2104 (51%)	108 (5.13%)	250 (11.88%)	-	-	1746 (82.98%)	-	1999 (48%)	104 (5.20%)	416 (20.81%)	-	-	1479 (73.99%)	-
Chow (US CDC)	6637	5143 (77%)	61 (1.19%)	80 (1.56%)	-	-	-	5002 (97.26%)	1494 (22%)	27 (1.81%)	78 (5.22%)	-	-	-	1389 (92.97%)
Argenziano	1000	151 (15%)	14 (9.27%)	18 (11.92%)	-	119 (78.81%)	-	-	849 (84%)	35 (4.12%)	161 (18.96%)	-	653 (76.91%)	-	-
Miyara	470	139 (29%)	14 (10.07%)	41 (29.50%)	-	77 (55.40%)	-	7 (5.04%)	340 (72%)	18 (5.29%)	111 (32.65%)	-	209 (61.47%)	-	2 (0.59%)
Lubetzky	54	15 (27%)	-	-	4 (26.67%)	-	-	11 (73.33%)	39 (72%)	-	-	8 (20.51%)	-	-	31 (79.49%)
Carrillo-Vega	9946	3922 (39%)	408 (10.40%)	-	-	-	-	3514 (89.60%)	6024 (60%)	486 (8.07%)	-	-	-	-	5538 (91.93%)
Yanover	4353	4180 (96%)	484 (11.58%)	118 (2.82%)	-	3578 (85.60%)	-	-	173 (3%)	30 (17.34%)	11 (6.36%)	-	132 (76.30%)	-	-
Hamer	387109	386349 (99%)	37333 (9.66%)	134542 (34.82%)	-	214474 (55.51%)	-	-	760 (0%)	93 (12.24%)	313 (41.18%)	-	354 (46.58%)	-	-
Heili-Frades	4712	1973 (41%)	121 (6.13%)	222 (11.25%)	-	-	1630 (82.62%)	1630 (82.62%)	2739 (58%)	112 (4.09%)	598 (21.83%)	-	-	2029 (74.08%)	-
Freites	123	69 (56%)	1 (1.45%)	-	-	-	-	68 (98.55%)	54 (43%)	3 (5.56%)	-	-	-	-	51 (94.44%)
Berumen	102875	18832 (18%)	-	-	1546 (8.21%)	-	17286 (91.79%)	-	12690 (12%)	-	-	1202 (9.47%)	-	11488 (90.53%)	-
Gianfrancesco	600	323 (53%)	-	-	61 (18.89%)	-	-	262 (81.11%)	277 (46%)	-	-	68 (24.55%)	-	-	209 (75.45%)
Chaudhry	40	19 (47%)	-	-	0 (0.00%)	-	-	19 (100.00%)	21 (52%)	-	-	6 (28.57%)	-	-	15 (71.43%)

nrow(table_3_word)

## [1] 14

save_as_docx(a, path = here('data_clean', 'Table_3.docx'))

0.5.3.1 Smoking and hospitalisation meta-analysis

The studies included in meta-analysis are c(“rentsch”, “argenziano”, “miyara_updated”, “yanover”, “hamer”).

# Data --------------------------------------------------------------------
meta <- tibble('author' = table_3$lead_author,
               'community_smoker' = table_3$community_current_smoker,
               'community_never_smoker' = table_3$community_never_smoker, 
               'hospitalised_smoker' = table_3$hospitalised_current_smoker, 
               'hospitalised_never_smoker' = table_3$hospitalised_never_smoker,
               'community_former_smoker' = table_3$community_former_smoker,
               'hospitalised_former_smoker' = table_3$hospitalised_former_smoker) %>%
  filter(author %in% included_studies$lead_author)
  
meta$author <- recode(meta$author, 'rentsch' = 'Rentsch',
                      'argenziano' = 'Argenziano',
                      'miyara_updated' = 'Miyara',
                      'yanover' = 'Yanover',
                      'hamer' = 'Hamer',
                      'heili_frades' = 'Heili-Frades')

# Current smoker hospitalisation ------------------------------------------
event_rates_smoker <- meta %>%
  mutate(., Ee = hospitalised_smoker) %>%
  mutate(., Ne = (hospitalised_smoker+community_smoker)) %>%
  mutate(., Ec = hospitalised_never_smoker) %>%
  mutate(., Nc = (hospitalised_never_smoker+community_never_smoker)) %>%
  rename('Author' = author) %>%
  select(Author, Ee, Ne, Ec, Nc)

a <- event_rates_smoker

event_rates_smoker <- metabin(Ee,
                              Ne,
                              Ec,
                              Nc,
                              data = event_rates_smoker,
                              studlab = paste(Author),
                              comb.fixed = F,
                              comb.random = T,
                              method.tau = 'SJ',
                              hakn = T,
                              prediction = F,
                              incr = 0.1,
                              sm = 'RR')
png(here::here('reports', 'figure', 'fig_5.png'), width=1480, height=546, res=120)
current_smoker_hospitalisation <- forest(event_rates_smoker,
                                sortvar = Author,
                                xlim = c(0.5, 3),
                                rightlabs = c('RR', '95% CI', 'Weight'),
                                leftlabs = c('Author', 'logRR', 'SE'),
                                lab.e = 'Current smoker',
                                lab.c = 'Never smoker',
                                print.tau2 = F,
                                col.diamond = 'blue',
                                col.diamond.lines = 'black',
                                col.square = 'black',
                                col.square.lines = 'black',
                                digits.sd = 2)
null <- dev.off()

a

## # A tibble: 5 x 5
##   Author        Ee    Ne    Ec     Nc
##   <chr>      <dbl> <dbl> <dbl>  <dbl>
## 1 Rentsch       90   159   106    216
## 2 Argenziano    35    49   653    772
## 3 Miyara        18    32   209    286
## 4 Yanover       30   514   132   3710
## 5 Hamer         93 37426   354 214828

Current versus never smokers

# Former smoker hospitalisation -------------------------------------------
event_rates_former <- meta %>%
  mutate(., Ee = hospitalised_former_smoker) %>%
  mutate(., Ne = (hospitalised_former_smoker+community_former_smoker)) %>%
  mutate(., Ec = hospitalised_never_smoker) %>%
  mutate(., Nc = (hospitalised_never_smoker+community_never_smoker)) %>%
  rename('Author' = author) %>%
  select(Author, Ee, Ne, Ec, Nc)

event_rates_former <- metabin(Ee,
                              Ne,
                              Ec,
                              Nc,
                              data = event_rates_former,
                              studlab = paste(Author),
                              comb.fixed = F,
                              comb.random = T,
                              method.tau = 'SJ',
                              hakn = T,
                              prediction = F,
                              incr = 0.1,
                              sm = 'RR')

a <- event_rates_former

png(here::here('reports', 'figure', 'fig_6.png'), width=1480, height=546, res=120)
former_smoker_hospitalisation <- forest(event_rates_former,
                                         sortvar = Author,
                                         xlim = c(0.5, 5),
                                         rightlabs = c('RR', '95% CI', 'Weight'),
                                         leftlabs = c('Author', 'logRR', 'SE'),
                                         lab.e = 'Former smoker',
                                         lab.c = 'Never smoker',
                                         print.tau2 = F,
                                         col.diamond = 'blue',
                                         col.diamond.lines = 'black',
                                         col.square = 'black',
                                         col.square.lines = 'black',
                                         digits.sd = 2)
null <- dev.off()
a

##                RR           95%-CI %W(random)
## Rentsch    1.0132 [0.8292; 1.2380]       20.9
## Argenziano 1.0634 [1.0039; 1.1263]       23.1
## Miyara     0.9993 [0.8867; 1.1262]       22.4
## Yanover    2.3966 [1.3292; 4.3214]       11.6
## Hamer      1.4085 [1.2100; 1.6396]       21.9
## 
## Number of studies combined: k = 5
## 
##                          RR           95%-CI    t p-value
## Random effects model 1.2134 [0.8234; 1.7883] 1.38  0.2383
## 
## Quantifying heterogeneity:
##  tau^2 = 0.0897 [0.0089; 1.0766]; tau = 0.2996 [0.0943; 1.0376];
##  I^2 = 81.0% [55.6%; 91.8%]; H = 2.29 [1.50; 3.50]
## 
## Test of heterogeneity:
##      Q d.f. p-value
##  21.02    4  0.0003
## 
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Sidik-Jonkman estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model

Former versus never smokers

0.5.4 Table 4: Smoking and COVID-19 disease severity

#Table 4
table_4_word <- table_4 %>%
  mutate(., sample = sample_with_severity) %>%
  mutate(., non_severe_disease_percentage = formatC(non_severe_disease/sample*100, digits = 2, format = "f")) %>%
  mutate(., non_severe_current_smoker_percent = formatC(non_severe_current_smoker/non_severe_disease*100, digits = 2, format = "f")) %>%
  mutate(., non_severe_former_smoker_percent = formatC(non_severe_former_smoker/non_severe_disease*100, digits = 2, format = "f")) %>%
  mutate(., non_severe_current_former_smoker_percent = formatC(non_severe_current_former_smoker/non_severe_disease*100, digits = 2, format = "f")) %>%
  mutate(., non_severe_never_smoker_percent = formatC(non_severe_never_smoker/non_severe_disease*100, digits = 2, format = "f")) %>%
  mutate(., non_severe_never_unknown_smoker_percent = formatC(non_severe_never_unknown_smoker/non_severe_disease*100, digits = 2, format = "f")) %>%
  mutate(., non_severe_not_stated_percent = formatC(non_severe_not_stated/non_severe_disease*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_number_percent = formatC(severe_disease_number/sample*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_current_smoker_percent = formatC(severe_disease_current_smoker/severe_disease_number*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_former_smoker_percent = formatC(severe_disease_former_smoker/severe_disease_number*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_current_former_smoker_percent = formatC(severe_disease_current_former_smoker/severe_disease_number*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_never_smoker_percent = formatC(severe_disease_never_smoker/severe_disease_number*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_never_unknown_percent = formatC(severe_disease_never_unknown/severe_disease_number*100, digits = 2, format = "f")) %>%
  mutate(., severe_disease_not_stated_percent = formatC(severe_disease_not_stated/severe_disease_number*100, digits = 2, format = "f")) %>%
  select(lead_author, sample, non_severe_disease, non_severe_disease_percentage, non_severe_current_smoker,
         non_severe_current_smoker_percent, non_severe_former_smoker, non_severe_former_smoker_percent,
         non_severe_current_former_smoker, non_severe_current_former_smoker_percent, non_severe_never_smoker,
         non_severe_never_smoker_percent, non_severe_never_unknown_smoker, non_severe_never_unknown_smoker_percent,
         non_severe_not_stated, non_severe_not_stated_percent, severe_disease_number, severe_disease_number_percent,
         severe_disease_current_smoker, severe_disease_current_smoker_percent, 
         severe_disease_former_smoker, severe_disease_former_smoker_percent, severe_disease_current_former_smoker,
         severe_disease_current_former_smoker_percent, severe_disease_never_smoker, severe_disease_never_smoker_percent,
         severe_disease_never_unknown, severe_disease_never_unknown_percent, severe_disease_not_stated,
         severe_disease_not_stated_percent) %>%
  write_rds(here::here('data_clean', 'table_4_word.rds'))

quality_table_4 <- table_4_word %>%
  left_join(., quality_rating, by = 'lead_author') %>%
  select(lead_author, overall_rating, missingness) %>%
  mutate(overall_quality = ifelse(overall_rating != 'poor' & missingness < 20, 'fair', 'poor')) %>%
  select(lead_author, overall_quality)

a <- table_4_word %>%
  filter(non_severe_disease >= 1) %>%
  mutate(Author = lead_author,
         Population = sample,
         non_severe_disease = paste(
           paste(non_severe_disease, 
                 as.integer(non_severe_disease_percentage), sep = " ("), "%)", sep = ""),
         ns_current_smoker = paste(paste(non_severe_current_smoker, (non_severe_current_smoker_percent), sep = " (")
                                  , "%)", sep = ""),
         ns_current_smoker = na_if(ns_current_smoker, "NA ( NA%)"),
         ns_former_smoker = paste(paste(non_severe_former_smoker, (non_severe_former_smoker_percent), sep = " (")
                                 , "%)", sep = ""),
         ns_former_smoker = na_if(ns_former_smoker, "NA ( NA%)"),
         ns_current_former_smoker = paste(paste(non_severe_current_former_smoker,
                                               (non_severe_current_former_smoker_percent), sep = " (")
                                         , "%)", sep = ""),
         ns_current_former_smoker = na_if(ns_current_former_smoker, "NA ( NA%)"),
         ns_never_smoker = paste(paste(non_severe_never_smoker,
                                      (non_severe_never_smoker_percent), sep = " (")
                                , "%)", sep = ""),
         ns_never_smoker = na_if(ns_never_smoker, "NA ( NA%)"),
         ns_never_unknown_smoker = paste(paste(non_severe_never_unknown_smoker,
                                              (non_severe_never_unknown_smoker_percent), sep = " (")
                                        , "%)", sep = ""),
         ns_never_unknown_smoker = na_if(ns_never_unknown_smoker, "NA ( NA%)"),
         ns_not_stated = paste(paste(non_severe_not_stated,
                                    (non_severe_not_stated_percent), sep = " (")
                              , "%)", sep = ""),
         ns_not_stated = na_if(ns_not_stated, "NA ( NA%)"),
         ) %>%
  mutate(severe_disease = paste(
    paste(severe_disease_number, 
          as.integer(severe_disease_number_percent), sep = " ("), "%)", sep = ""),
    s_current_smoker = paste(paste(severe_disease_current_smoker, (severe_disease_current_smoker_percent), sep = " (")
                             , "%)", sep = ""),
    s_current_smoker = na_if(s_current_smoker, "NA ( NA%)"),
    s_former_smoker = paste(paste(severe_disease_former_smoker, (severe_disease_former_smoker_percent), sep = " (")
                            , "%)", sep = ""),
    s_former_smoker = na_if(s_former_smoker, "NA ( NA%)"),
    s_current_former_smoker = paste(paste(severe_disease_current_former_smoker,
                                          (severe_disease_current_former_smoker_percent), sep = " (")
                                    , "%)", sep = ""),
    s_current_former_smoker = na_if(s_current_former_smoker, "NA ( NA%)"),
    s_never_smoker = paste(paste(severe_disease_never_smoker,
                                 (severe_disease_never_smoker_percent), sep = " (")
                           , "%)", sep = ""),
    s_never_smoker = na_if(s_never_smoker, "NA ( NA%)"),
    s_never_unknown_smoker = paste(paste(severe_disease_never_unknown,
                                         (severe_disease_never_unknown_percent), sep = " (")
                                   , "%)", sep = ""),
    s_never_unknown_smoker = na_if(s_never_unknown_smoker, "NA ( NA%)"),
    s_not_stated = paste(paste(severe_disease_not_stated,
                                 (severe_disease_not_stated_percent), sep = " (")
                           , "%)", sep = ""),
    s_not_stated = na_if(s_not_stated, "NA ( NA%)")) %>%
  select(Author, Population, non_severe_disease, ns_current_smoker, ns_former_smoker,
         ns_current_former_smoker, ns_never_smoker, ns_never_unknown_smoker, ns_not_stated,
         severe_disease, s_current_smoker, s_former_smoker, s_current_former_smoker, s_never_smoker, s_never_unknown_smoker,
         s_not_stated)

         
a$Author <- to_upper_camel_case(a$Author, sep_out = ", ")
a$Author  <- a$Author  %>%
  recode("De, La, Rica"  = "de la Rica",
         "Mejia, Vilet" = "Mejia-Vilet",
         "Vaquero, Roncero" = "Vaquero-Roncero")

numeric_columns <- 'Population'

a <- flextable(a) %>%
    set_caption(caption = 'COVID-19 severity by smoking status') %>%
  colformat_num(col_keys = numeric_columns, digits = 0, na_str = '-', big.mark = ',') %>%
  colformat_char(na_str = '-')

a <- set_header_labels(a,
                       Population = "Population with severity",
                       non_severe_disease = "N (%)",
                       ns_current_smoker = "Current smoker (%)",
                       ns_former_smoker = "Former smoker (%)",
                       ns_current_former_smoker = "Current/former smoker (%)",
                       ns_never_unknown_smoker = "Never/unknown smoker (%)",
                       ns_never_smoker = "Never smoker (%)",
                       ns_not_stated = "Not stated (%)",
                       severe_disease = "N (%)", 
                       s_current_smoker = "Current smoker (%)", 
                       s_former_smoker = "Former smoker (%)", 
                       s_current_former_smoker = "Current/former smoker (%)",
                       s_never_unknown_smoker = "Never/unknown smoker (%)",
                       s_never_smoker = "Never smoker (%)",
                       s_not_stated = "Not stated (%)") %>%
  add_header_row(top = TRUE, values = c("","Non severe disease", "Severe disease" ), colwidths = c(2, 7, 7)) %>%
  theme_booktabs() %>%
  fix_border_issues() %>%
  set_table_properties(width = 1, layout = 'autofit')
a

COVID-19 severity by smoking status

		Non severe disease							Severe disease
Author	Population with severity	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Never/unknown smoker (%)	Not stated (%)	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Never/unknown smoker (%)	Not stated (%)
Guan, Ni	1085	913 (84%)	108 (11.83%)	12 (1.31%)	-	793 (86.86%)	-	-	172 (15%)	29 (16.86%)	9 (5.23%)	-	134 (77.91%)	-	-
Zhang, Dong	9	3 (33%)	0 (0.00%)	3 (100.00%)	-	0 (0.00%)	-	-	6 (66%)	2 (33.33%)	4 (66.67%)	-	0 (0.00%)	-	-
Wan	9	8 (88%)	8 (100.00%)	0 (0.00%)	-	0 (0.00%)	-	-	1 (11%)	1 (100.00%)	0 (0.00%)	-	0 (0.00%)	-	-
Huang, Wang	3	3 (100%)	3 (100.00%)	0 (0.00%)	-	0 (0.00%)	-	-	0 (0%)	0 (NaN%)	0 (NaN%)	-	0 (NaN%)	-	-
Rentsch	285	168 (58%)	47 (27.98%)	53 (31.55%)	-	68 (40.48%)	-	-	117 (41%)	43 (36.75%)	36 (30.77%)	-	38 (32.48%)	-	-
Hu	323	151 (46%)	-	-	12 (7.95%)	-	139 (92.05%)	-	172 (53%)	-	-	26 (15.12%)	-	146 (84.88%)	-
Wang, Pan	125	100 (80%)	-	-	9 (9.00%)	-	91 (91.00%)	-	25 (20%)	-	-	7 (28.00%)	-	18 (72.00%)	-
Petrilli	1278	932 (72%)	62 (6.65%)	175 (18.78%)	-	-	391 (41.95%)	-	650 (50%)	28 (4.31%)	145 (22.31%)	-	-	477 (73.38%)	-
Kim	27	21 (77%)	3 (14.29%)	-	-	-	18 (85.71%)	-	6 (22%)	2 (33.33%)	0 (0.00%)	-	-	4 (66.67%)	-
Shi, Yu	474	425 (89%)	-	-	34 (8.00%)	-	391 (92.00%)	-	49 (10%)	-	-	6 (12.24%)	-	43 (87.76%)	-
Liao, Feng	148	92 (62%)	-	-	5 (5.43%)	-	-	87 (94.57%)	56 (37%)	3 (5.36%)	-	-	-	-	53 (94.64%)
Shi, Ren	134	88 (65%)	-	-	8 (9.09%)	-	-	80 (90.91%)	46 (34%)	-	-	6 (13.04%)	-	-	40 (86.96%)
Hadjadj	50	15 (30%)	1 (6.67%)	2 (13.33%)	-	12 (80.00%)	-	-	35 (70%)	0 (0.00%)	7 (20.00%)	-	28 (80.00%)	-	-
Zheng, Xiong	73	43 (58%)	-	-	6 (13.95%)	37 (86.05%)	-	-	30 (41%)	-	-	2 (6.67%)	28 (93.33%)	-	-
de la Rica	48	26 (54%)	-	-	6 (23.08%)	-	-	20 (76.92%)	20 (41%)	-	-	4 (20.00%)	-	-	16 (80.00%)
Yin, Yang	106	47 (44%)	-	-	6 (12.77%)	-	-	41 (87.23%)	59 (55%)	-	-	12 (20.34%)	-	-	47 (79.66%)
Allenbach	147	100 (68%)	-	-	9 (9.00%)	-	-	91 (91.00%)	47 (31%)	-	-	0 (0.00%)	-	-	47 (100.00%)
Goyal	393	263 (66%)	14 (5.32%)	-	-	-	-	249 (94.68%)	130 (33%)	6 (4.62%)	-	-	-	-	124 (95.38%)
Feng	454	333 (73%)	27 (8.11%)	-	-	-	-	306 (91.89%)	121 (26%)	17 (14.05%)	-	-	-	-	104 (85.95%)
Yao	108	83 (76%)	1 (1.20%)	-	-	-	-	82 (98.80%)	25 (23%)	3 (12.00%)	-	-	-	-	22 (88.00%)
Sami	490	400 (81%)	53 (13.25%)	-	-	-	-	347 (86.75%)	90 (18%)	16 (17.78%)	-	-	-	-	74 (82.22%)
Regina	200	163 (81%)	9 (5.52%)	-	-	-	-	154 (94.48%)	37 (18%)	0 (0.00%)	-	-	-	-	37 (100.00%)
Feuth	28	21 (75%)	1 (4.76%)	7 (33.33%)	-	13 (61.90%)	-	-	7 (25%)	2 (28.57%)	1 (14.29%)	-	4 (57.14%)	-	-
Mejia-Vilet	329	214 (65%)	-	-	13 (6.07%)	-	-	201 (93.93%)	115 (34%)	-	-	10 (8.70%)	-	-	105 (91.30%)
Chen, Jiang	135	54 (40%)	-	-	4 (7.41%)	-	-	50 (92.59%)	81 (60%)	-	-	9 (11.11%)	-	-	72 (88.89%)
Vaquero-Roncero	146	75 (51%)	-	-	4 (5.33%)	-	-	71 (94.67%)	71 (48%)	-	-	6 (8.45%)	-	-	65 (91.55%)
Kim, Garg	2490	1692 (67%)	112 (6.62%)	395 (23.35%)	-	-	1185 (70.04%)	-	798 (32%)	38 (4.76%)	247 (30.95%)	-	-	512 (64.16%)	-
Wu	174	92 (52%)	-	-	47 (51.09%)	-	45 (48.91%)	-	82 (47%)	11 (13.41%)	-	-	-	71 (86.59%)	-
Russell	106	87 (82%)	6 (6.90%)	22 (25.29%)	-	35 (40.23%)	-	24 (27.59%)	19 (17%)	1 (5.26%)	4 (21.05%)	-	7 (36.84%)	-	7 (36.84%)
Chaudhry	40	34 (85%)	-	-	5 (14.71%)	-	-	29 (85.29%)	6 (15%)	-	-	1 (16.67%)	-	-	5 (83.33%)
Patel	104	67 (64%)	25 (37.31%)	-	-	-	36 (53.73%)	6 (8.96%)	37 (35%)	18 (48.65%)	-	-	-	15 (40.54%)	4 (10.81%)
Garibaldi	832	532 (63%)	25 (4.70%)	107 (20.11%)	-	-	-	400 (75.19%)	300 (36%)	21 (7.00%)	81 (27.00%)	-	-	-	198 (66.00%)

nrow(table_4_word)

## [1] 33

save_as_docx(a, path = here('data_clean', 'Table_4.docx'))

0.5.4.1 Smoking and severity meta-analysis

The studies included in meta-analysis are c(“guan_ni”, “rentsch”, “hadjadj”, “feuth”).

# Data --------------------------------------------------------------------
meta <- tibble('author' = table_4$lead_author,
               'non_severe_smoker' = table_4$non_severe_current_smoker,
               'non_severe_never_smoker' = table_4$non_severe_never_smoker,
               'severe_smoker' = table_4$severe_disease_current_smoker,
               'severe_never_smoker' = table_4$severe_disease_never_smoker,
               'non_severe_former_smoker' = table_4$non_severe_former_smoker,
               'severe_former_smoker' = table_4$severe_disease_former_smoker) %>%
  filter(author %in% included_studies$lead_author)

meta$author <- recode(meta$author, 'rentsch' = 'Rentsch',
                      'guan_ni' = 'Guan, Ni',
                      'hadjadj' = 'Hadjadj',
                      "feuth" = "Feuth")


# Current smoker severity ------------------------------------------
event_rates_smoker <- meta %>%
  mutate(., Ee = severe_smoker) %>%
  mutate(., Ne = (severe_smoker+non_severe_smoker)) %>%
  mutate(., Ec = severe_never_smoker) %>%
  mutate(., Nc = (severe_never_smoker+non_severe_never_smoker)) %>%
  rename('Author' = author) %>%
  select(Author, Ee, Ne, Ec, Nc)

event_rates_smoker <- metabin(Ee,
                              Ne,
                              Ec,
                              Nc,
                              data = event_rates_smoker,
                              studlab = paste(Author),
                              comb.fixed = T,
                              comb.random = F,
                              method.tau = 'SJ',
                              hakn = F,
                              prediction = F,
                              incr = 0.1,
                              sm = 'RR')
png(here::here('reports', 'figure', 'fig_7.png'), width=1024, height=546, res=120)
current_smoker_severity <- forest(event_rates_smoker,
                                         sortvar = Author,
                                         xlim = c(0.1, 5),
                                         rightlabs = c('RR', '95% CI', 'Weight'),
                                         leftlabs = c('Author', 'logRR', 'SE'),
                                         lab.e = 'Current smoker',
                                         lab.c = 'Never smoker',
                                         print.tau2 = F,
                                         col.diamond = 'blue',
                                         col.diamond.lines = 'black',
                                         col.square = 'black',
                                         col.square.lines = 'black',
                                         digits.sd = 2)
null <- dev.off()

event_rates_smoker

##              RR            95%-CI %W(fixed)
## Guan, Ni 1.4644 [1.0226;  2.0970]      47.8
## Rentsch  1.3327 [0.9544;  1.8611]      48.3
## Hadjadj  0.1297 [0.0004; 47.9753]       2.3
## Feuth    2.8333 [0.8772;  9.1514]       1.7
## 
## Number of studies combined: k = 4
## 
##                        RR           95%-CI    z p-value
## Fixed effect model 1.3932 [1.0947; 1.7731] 2.70  0.0070
## 
## Quantifying heterogeneity:
##  tau^2 = 0.3716 [0.0000; 18.9140]; tau = 0.6096 [0.0000; 4.3490];
##  I^2 = 0.0% [0.0%; 78.4%]; H = 1.00 [1.00; 2.15]
## 
## Test of heterogeneity:
##     Q d.f. p-value
##  2.13    3  0.5465
## 
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Sidik-Jonkman estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Continuity correction of 0.1 in studies with zero cell frequencies

Current versus never smokers

# Former smoker hospitalisation -------------------------------------------
event_rates_former <- meta %>%
  mutate(., Ee = severe_former_smoker) %>%
  mutate(., Ne = (severe_former_smoker+non_severe_former_smoker)) %>%
  mutate(., Ec = severe_never_smoker) %>%
  mutate(., Nc = (severe_never_smoker+non_severe_never_smoker)) %>%
  rename('Author' = author) %>%
  select(Author, Ee, Ne, Ec, Nc)

event_rates_former <- metabin(Ee,
                              Ne,
                              Ec,
                              Nc,
                              data = event_rates_former,
                              studlab = paste(Author),
                              comb.fixed = F,
                              comb.random = T,
                              method.tau = 'SJ',
                              hakn = F,
                              prediction = F,
                              incr = 0.1,
                              sm = 'RR')

png(here::here('reports', 'figure', 'fig_8.png'), width=1024, height=546, res=120)
former_smoker_severity <- forest(event_rates_former,
                                        sortvar = Author,
                                        xlim = c(0.2, 10),
                                        rightlabs = c('RR', '95% CI', 'Weight'),
                                        leftlabs = c('Author', 'logRR', 'SE'),
                                        lab.e = 'Former smoker',
                                        lab.c = 'Never smoker',
                                        print.tau2 = F,
                                        col.diamond = 'blue',
                                        col.diamond.lines = 'black',
                                        col.square = 'black',
                                        col.square.lines = 'black',
                                        digits.sd = 2)
null <- dev.off()

event_rates_former

##              RR           95%-CI %W(random)
## Guan, Ni 2.9648 [1.7660; 4.9774]       28.9
## Rentsch  1.1283 [0.7885; 1.6146]       32.4
## Hadjadj  1.1111 [0.7419; 1.6640]       31.4
## Feuth    0.5312 [0.0702; 4.0199]        7.3
## 
## Number of studies combined: k = 4
## 
##                          RR           95%-CI    z p-value
## Random effects model 1.4045 [0.7605; 2.5938] 1.09  0.2778
## 
## Quantifying heterogeneity:
##  tau^2 = 0.2693 [0.0153; 6.4303]; tau = 0.5189 [0.1238; 2.5358];
##  I^2 = 74.0% [27.2%; 90.7%]; H = 1.96 [1.17; 3.28]
## 
## Test of heterogeneity:
##      Q d.f. p-value
##  11.53    3  0.0092
## 
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Sidik-Jonkman estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau

Former versus never smokers

0.5.5 Table 5: Smoking and COVID-19 mortality

#Table 5
table_5_word <- table_5 %>%
  mutate(., sample = sample_with_deaths) %>%
  mutate(., deaths_percentage = formatC(deaths/sample*100, digits = 2, format = "f")) %>%
  mutate(., death_current_smokers_percent = formatC(death_current_smokers/deaths*100, digits = 2, format = "f")) %>%
  mutate(., death_former_smokers_percent = formatC(death_former_smokers/deaths*100, digits = 2, format = "f")) %>%
  mutate(., death_current_former_smokers_percent = formatC(death_current_former_smokers/deaths*100, digits = 2, format = "f")) %>%
  mutate(., death_never_smokers_percent = formatC(death_never_smokers/deaths*100, digits = 2, format = "f")) %>%
  mutate(., death_never_unknown_smokers_percent = formatC(death_never_unknown_smokers/deaths*100, digits = 2, format = "f")) %>%
  mutate(., death_not_stated_percent = formatC(death_not_stated/deaths*100, digits = 2, format = "f")) %>%
  mutate(., recovered_percentage = formatC(recovered/sample*100, digits = 2, format = "f")) %>%
  mutate(., recovered_current_smokers_percent = formatC(recovered_current_smoking/recovered*100, digits = 2, format = "f")) %>%
  mutate(., recovered_former_smokers_percent = formatC(recovered_former_smoker/recovered*100, digits = 2, format = "f")) %>%
  mutate(., recovered_current_former_smokers_percent = formatC(recovered_current_former_smokers/recovered*100, digits = 2, format = "f")) %>%
  mutate(., recovered_never_smokers_percent = formatC(recovered_never_smoker/recovered*100, digits = 2, format = "f")) %>%
  mutate(., recovered_never_unknown_smokers_percent = formatC(recovered_never_unknown_smoker/recovered*100, digits = 2, format = "f")) %>%
  mutate(., recovered_not_stated_percent = formatC(recovered_not_stated/recovered*100, digits = 2, format = "f")) %>%
  select(lead_author, sample, recovered, recovered_percentage, recovered_current_smoking, recovered_current_smokers_percent,
         recovered_former_smoker, recovered_former_smokers_percent, recovered_current_former_smokers,
         recovered_current_former_smokers_percent, recovered_never_smoker, recovered_never_smokers_percent,
         recovered_never_unknown_smoker, recovered_never_unknown_smokers_percent, recovered_not_stated,
         recovered_not_stated_percent,
         deaths, deaths_percentage, death_current_smokers, death_current_smokers_percent,
         death_former_smokers, death_former_smokers_percent, death_current_former_smokers,
         death_current_former_smokers_percent, death_never_smokers, death_never_smokers_percent,
         death_never_unknown_smokers, death_never_unknown_smokers_percent, death_not_stated, death_not_stated_percent) %>%
  write_rds(here::here('data_clean', 'table_5_word.rds'))

quality_table_5 <- table_5_word %>%
  left_join(., quality_rating, by = 'lead_author') %>%
  select(lead_author, overall_rating, missingness) %>%
  mutate(overall_quality = ifelse(overall_rating != 'poor' & missingness < 20, 'fair', 'poor')) %>%
  select(lead_author, overall_quality)

a <- table_5_word %>%
  filter(recovered >= 1) %>%
  mutate(Author = lead_author,
         Population = sample,
         recovered = paste(
           paste(recovered, 
                 as.integer(recovered_percentage), sep = " ("), "%)", sep = ""),
         recovered_current_smoker = paste(paste(recovered_current_smoking, 
                                                 (recovered_current_smokers_percent), sep = " (")
                                  , "%)", sep = ""),
         recovered_current_smoker = na_if(recovered_current_smoker, "NA ( NA%)"),
         recovered_former_smoker = paste(paste(recovered_former_smoker,
                                               (recovered_former_smokers_percent), sep = " (")
                                 , "%)", sep = ""),
         recovered_former_smoker = na_if(recovered_former_smoker, "NA ( NA%)"),
         recovered_current_former_smoker = paste(paste(recovered_current_former_smokers,
                                               (recovered_current_former_smokers_percent), sep = " (")
                                         , "%)", sep = ""),
         recovered_current_former_smoker = na_if(recovered_current_former_smoker, "NA ( NA%)"),
         recovered_never_smoker = paste(paste(recovered_never_smoker,
                                      (recovered_never_smokers_percent), sep = " (")
                                , "%)", sep = ""),
         recovered_never_smoker = na_if(recovered_never_smoker, "NA ( NA%)"),
         recovered_never_unknown_smoker = paste(paste(recovered_never_unknown_smoker,
                                              (recovered_never_unknown_smokers_percent), sep = " (")
                                        , "%)", sep = ""),
         recovered_never_unknown_smoker = na_if(recovered_never_unknown_smoker, "NA ( NA%)"),
         recovered_not_stated = paste(paste(recovered_not_stated,
                                    (recovered_not_stated_percent), sep = " (")
                              , "%)", sep = ""),
         recovered_not_stated = na_if(recovered_not_stated, "NA ( NA%)")) %>%
  mutate(deaths = paste(
           paste(deaths, 
                 as.integer(deaths_percentage), sep = " ("), "%)", sep = ""),
         death_current_smoker = paste(paste(death_current_smokers, 
                                                 (death_current_smokers_percent), sep = " (")
                                  , "%)", sep = ""),
         death_current_smoker = na_if(death_current_smoker, "NA ( NA%)"),
         death_former_smoker = paste(paste(death_former_smokers,
                                               (death_former_smokers_percent), sep = " (")
                                 , "%)", sep = ""),
         death_former_smoker = na_if(death_former_smoker, "NA ( NA%)"),
         death_current_former_smoker = paste(paste(death_current_former_smokers,
                                               (death_current_former_smokers_percent), sep = " (")
                                         , "%)", sep = ""),
         death_current_former_smoker = na_if(death_current_former_smoker, "NA ( NA%)"),
         death_never_smoker = paste(paste(death_never_smokers,
                                      (death_never_smokers_percent), sep = " (")
                                , "%)", sep = ""),
         death_never_smoker = na_if(death_never_smoker, "NA ( NA%)"),
         death_never_unknown_smoker = paste(paste(death_never_unknown_smokers,
                                              (death_never_unknown_smokers_percent), sep = " (")
                                        , "%)", sep = ""),
         death_never_unknown_smoker = na_if(death_never_unknown_smoker, "NA ( NA%)"),
         death_not_stated = paste(paste(death_not_stated,
                                    (death_not_stated_percent), sep = " (")
                              , "%)", sep = ""),
         death_not_stated = na_if(death_not_stated, "NA ( NA%)")) %>%
  select(Author, Population, recovered, recovered_current_smoker, recovered_former_smoker, recovered_current_former_smoker,
         recovered_never_smoker, recovered_never_unknown_smoker, recovered_not_stated, deaths, death_current_smoker,
         death_former_smoker, death_current_former_smoker, death_never_smoker, death_never_unknown_smoker,
         death_not_stated)

         
a$Author <- to_upper_camel_case(a$Author, sep_out = ", ")
a$Author  <- a$Author  %>%
  recode("Carrillo, Vega"  = "Carillo-Vega",
         "Heili, Frades" = "Heili-Frades",
         "Al, Hindawi" = "Al-Hindawi",
         "Soto, Mota" = "Soto-Mota")

numeric_columns <- 'Population'

a <- flextable(a) %>%
    set_caption(caption = 'COVID-19 mortality by smoking status') %>%
  colformat_num(col_keys = numeric_columns, digits = 0, na_str = '-', big.mark = ',') %>%
  colformat_char(na_str = '-')

a <- set_header_labels(a,
                       Population = "Population with mortality",
                       recovered = "N (%)",
                       recovered_current_smoker = "Current smoker (%)",
                       recovered_former_smoker = "Former smoker (%)",
                       recovered_current_former_smoker = "Current/former smoker (%)",
                       recovered_never_unknown_smoker = "Never/unknown smoker (%)",
                       recovered_never_smoker = "Never smoker (%)",
                       recovered_not_stated = "Not stated (%)",
                       deaths = "N (%)", 
                       death_current_smoker = "Current smoker (%)", 
                       death_former_smoker = "Former smoker (%)", 
                       death_current_former_smoker = "Current/former smoker (%)",
                       death_never_unknown_smoker = "Never/unknown smoker (%)",
                       death_never_smoker = "Never smoker (%)",
                       death_not_stated = "Not stated (%)") %>%
  add_header_row(top = TRUE, values = c("","Recovered", "Died" ), colwidths = c(2, 7, 7)) %>%
  theme_booktabs() %>%
  fix_border_issues() %>%
  set_table_properties(width = 1, layout = 'autofit')
a

COVID-19 mortality by smoking status

		Recovered							Died
Author	Population with mortality	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Never/unknown smoker (%)	Not stated (%)	N (%)	Current smoker (%)	Former smoker (%)	Current/former smoker (%)	Never smoker (%)	Never/unknown smoker (%)	Not stated (%)
Chen	274	161 (58%)	5 (3.11%)	5 (3.11%)	-	-	-	151 (93.79%)	113 (41%)	7 (6.19%)	2 (1.77%)	-	-	-	104 (92.04%)
Zhou, Yu	191	137 (71%)	6 (4.38%)	-	-	-	-	131 (95.62%)	54 (28%)	5 (9.26%)	-	-	-	-	49 (90.74%)
Yang, Yu	52	20 (38%)	2 (10.00%)	-	-	-	18 (90.00%)	-	32 (61%)	-	-	-	-	32 (100.00%)	-
Gaibazzi	441	285 (64%)	11 (3.86%)	30 (10.53%)	-	244 (85.61%)	-	-	156 (35%)	10 (6.41%)	14 (8.97%)	-	132 (84.62%)	-	-
Borobia	2226	1766 (79%)	113 (6.40%)	-	-	-	-	1653 (93.60%)	460 (20%)	44 (9.57%)	-	-	-	-	416 (90.43%)
Giacomelli	233	185 (79%)	-	-	53 (28.65%)	132 (71.35%)	-	-	48 (20%)	-	-	17 (35.42%)	31 (64.58%)	-	0 (0.00%)
Yao	108	96 (88%)	1 (1.04%)	-	-	-	-	95 (98.96%)	12 (11%)	3 (25.00%)	-	-	-	-	9 (75.00%)
Carillo-Vega	9946	8983 (90%)	795 (8.85%)	-	-	-	-	8188 (91.15%)	963 (9%)	99 (10.28%)	-	-	-	-	864 (89.72%)
Heng	51	39 (76%)	6 (15.38%)	-	-	-	-	33 (84.62%)	12 (23%)	1 (8.33%)	-	-	-	-	11 (91.67%)
Heili-Frades	4712	4086 (86%)	210 (5.14%)	659 (16.13%)	-	-	3217 (78.73%)	-	626 (13%)	23 (3.67%)	161 (25.72%)	-	-	442 (70.61%)	-
Kim, Garg	2490	2070 (83%)	128 (6.18%)	481 (23.24%)	-	-	1461 (70.58%)	-	420 (16%)	22 (5.24%)	161 (38.33%)	-	-	236 (56.19%)	-
Al-Hindawi	31	15 (48%)	0 (0.00%)	10 (66.67%)	-	5 (33.33%)	-	-	16 (51%)	1 (6.25%)	12 (75.00%)	-	3 (18.75%)	-	-
Louis	22	16 (72%)	-	-	7 (43.75%)	-	-	9 (56.25%)	6 (27%)	-	-	3 (50.00%)	-	-	3 (50.00%)
Soto-Mota	400	200 (50%)	-	-	23 (11.50%)	-	-	177 (88.50%)	200 (50%)	-	-	25 (12.50%)	-	-	175 (87.50%)
Garibaldi	747	634 (84%)	36 (5.68%)	129 (20.35%)	-	-	-	469 (73.97%)	113 (15%)	6 (5.31%)	36 (31.86%)	-	-	-	71 (62.83%)
Docherty	13364	8199 (61%)	370 (4.51%)	1832 (22.34%)	-	4179 (50.97%)	-	1818 (22.17%)	5165 (38%)	214 (4.14%)	1350 (26.14%)	-	2105 (40.76%)	-	1496 (28.96%)

nrow(table_5_word)

## [1] 17

save_as_docx(a, path = here('data_clean', 'Table_5.docx'))

0.5.5.1 Smoking and COVID-19 mortality meta-analysis

The studies included in meta-analysis are c(“gaibazzi”, “al_hindawi”).

# Data --------------------------------------------------------------------
meta <- tibble('author' = table_5$lead_author,
               'recovered_current_smoker' = table_5$recovered_current_smoking,
               'recovered_never_smoker' = table_5$recovered_never_smoker,
               'death_current_smoker' = table_5$death_current_smokers,
               'death_never_smoker' = table_5$death_never_smokers,
               'recovered_former_smoker' = table_5$recovered_former_smoker,
               'death_former_smoker' = table_5$death_former_smokers) %>%
  filter(author %in% included_studies$lead_author)

meta$author <- recode(meta$author, 'gaibazzi' = 'Gaibazzi',
                      'al_hindawi' = 'Al-Hindawi')


# Current smoker mortality ------------------------------------------
event_rates_smoker <- meta %>%
  mutate(., Ee = death_current_smoker) %>%
  mutate(., Ne = (death_current_smoker+recovered_current_smoker)) %>%
  mutate(., Ec = death_never_smoker) %>%
  mutate(., Nc = (death_never_smoker+recovered_never_smoker)) %>%
  rename('Author' = author) %>%
  select(Author, Ee, Ne, Ec, Nc)

event_rates_smoker <- metabin(Ee,
                              Ne,
                              Ec,
                              Nc,
                              data = event_rates_smoker,
                              studlab = paste(Author),
                              comb.fixed = T,
                              comb.random = F,
                              method.tau = 'SJ',
                              hakn = F,
                              prediction = F,
                              incr = 0.1,
                              sm = 'RR')
png(here::here('reports', 'figure', 'fig_9.png'), width=1024, height=546, res=120)
current_smoker_mortality <- forest(event_rates_smoker,
                                         sortvar = Author,
                                         xlim = c(0.1, 5),
                                         rightlabs = c('RR', '95% CI', 'Weight'),
                                         leftlabs = c('Author', 'logRR', 'SE'),
                                         lab.e = 'Current smoker',
                                         lab.c = 'Never smoker',
                                         print.tau2 = F,
                                         col.diamond = 'blue',
                                         col.diamond.lines = 'black',
                                         col.square = 'black',
                                         col.square.lines = 'black',
                                         digits.sd = 2)
null <- dev.off()

event_rates_smoker

##                RR           95%-CI %W(fixed)
## Gaibazzi   1.3564 [0.8485; 2.1684]      94.6
## Al-Hindawi 2.6129 [1.0897; 6.2655]       5.4
## 
## Number of studies combined: k = 2
## 
##                        RR           95%-CI    z p-value
## Fixed effect model 1.4137 [0.9097; 2.1970] 1.54  0.1238
## 
## Quantifying heterogeneity:
##  tau^2 = 0.0980; tau = 0.3130; I^2 = 40.3%; H = 1.29
## 
## Test of heterogeneity:
##     Q d.f. p-value
##  1.68    1  0.1954
## 
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Sidik-Jonkman estimator for tau^2
## - Continuity correction of 0.1 in studies with zero cell frequencies

Current versus never smokers

# Former smoker mortality ------------------------------------------
event_rates_former <- meta %>%
  mutate(., Ee = death_former_smoker) %>%
  mutate(., Ne = (death_former_smoker+recovered_former_smoker)) %>%
  mutate(., Ec = death_never_smoker) %>%
  mutate(., Nc = (death_never_smoker+recovered_never_smoker)) %>%
  rename('Author' = author) %>%
  select(Author, Ee, Ne, Ec, Nc)

event_rates_former <- metabin(Ee,
                              Ne,
                              Ec,
                              Nc,
                              data = event_rates_former,
                              studlab = paste(Author),
                              comb.fixed = T,
                              comb.random = F,
                              method.tau = 'SJ',
                              hakn = F,
                              prediction = F,
                              incr = 0.1,
                              sm = 'RR')
png(here::here('reports', 'figure', 'fig_10.png'), width=1024, height=546, res=120)
former_smoker_mortality <- forest(event_rates_former,
                                         sortvar = Author,
                                         xlim = c(0.1, 5),
                                         rightlabs = c('RR', '95% CI', 'Weight'),
                                         leftlabs = c('Author', 'logRR', 'SE'),
                                         lab.e = 'Former smoker',
                                         lab.c = 'Never smoker',
                                         print.tau2 = F,
                                         col.diamond = 'blue',
                                         col.diamond.lines = 'black',
                                         col.square = 'black',
                                         col.square.lines = 'black',
                                         digits.sd = 2)
null <- dev.off()

event_rates_former

##                RR           95%-CI %W(fixed)
## Gaibazzi   0.9063 [0.5757; 1.4269]      86.3
## Al-Hindawi 1.4545 [0.5500; 3.8467]      13.7
## 
## Number of studies combined: k = 2
## 
##                        RR           95%-CI     z p-value
## Fixed effect model 0.9816 [0.6515; 1.4788] -0.09  0.9292
## 
## Quantifying heterogeneity:
##  tau^2 = 0.0304; tau = 0.1744; I^2 = 0.0%; H = 1.00
## 
## Test of heterogeneity:
##     Q d.f. p-value
##  0.75    1  0.3876
## 
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Sidik-Jonkman estimator for tau^2

Former versus never smokers

0.5.6 Study quality assessment

#Table 6
table_6_word <- table_6 %>%
  write_rds(here::here('data_clean', 'table_6_word.rds'))

a <- table_6_word %>%
  rename("Lead Author" = 1,
         "Question 1." = 2,
         "Question 2." = 3,
         "Question 3." = 4,
         "Question 4." = 5,
         "Question 5." = 6,
         "Question 6." = 7,
         "Question 7." = 8,
         "Question 8." = 9,
         "Question 9." = 10,
         "Question 10." = 11,
         "Question 11." = 12,
         "Question 12." = 13,
         "Question 13." = 14,
         "Question 14." = 15,
         "Missing data" = 17,
         "Overall rating (Poor, Fair, Good)" = 16) %>%
  mutate(`Missing data` = if_else(is.na(`Missing data`), 0, `Missing data`)) %>%
  select(1:15, 17, 16)

a$`Lead Author` <- to_upper_camel_case(a$`Lead Author`, sep_out = ", ")
a$`Lead Author`  <- a$`Lead Author`  %>%
  recode("Carrillo, Vega"  = "Carillo-Vega",
         "Heili, Frades" = "Heili-Frades",
         "Al, Hindawi" = "Al-Hindawi",
         "Soto, Mota" = "Soto-Mota",
         "Chow, Us, Cdc" = "Chow (US CDC)",
         "Gold, Us, Cdc" = "Gold (US CDC)",
         "Miyara, Updated" = "Miyars",
         "De, La, Rica" = "de la Rica",
         "Opensafely, Collaborative" = "The Opensafely Collaborative",
         "Bello, Chavolla" = "Bello-Chavolla",
         "De, Lusignan" = "de Lusignan",
         "Mejia, Vilet" = "Megia-Vilet",
         "Vaquero, Roncero" = "Vaquero-Roncer",
         "Del, Valle" = "del Valle")
a <- a %>%
  mutate("Missing data" = formatC(a$`Missing data`, digits = 2, format = "f"))
a$`Missing data` <- paste0(a$`Missing data`, "%")
         
         
         
a <- flextable(a) %>%
    set_caption(caption = 'National Institutes of Health’s Quality Assessment Tool for Observational Cohort and Cross-Sectional Studies') %>%
  colformat_char(na_str = '-') %>%
  theme_booktabs() %>%
  fix_border_issues() %>%
  set_table_properties(width = 1, layout = 'autofit')
a

National Institutes of Health’s Quality Assessment Tool for Observational Cohort and Cross-Sectional Studies

Lead Author	Question 1.	Question 2.	Question 3.	Question 4.	Question 5.	Question 6.	Question 7.	Question 8.	Question 9.	Question 10.	Question 11.	Question 12.	Question 13.	Question 14.	Missing data	Overall rating (Poor, Fair, Good)
Guan, Ni	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	1.27%	Fair
Guan, Liang	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	0.00%	Poor
Lian	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	No	No	Not applicable	No	93.15%	Poor
Jin	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.70%	Poor
Chen	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.07%	Poor
Zhou, Yu	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	94.24%	Poor
Mo	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	96.13%	Poor
Zhang, Dong	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.57%	Poor
Wan	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.33%	Poor
Liu, Tao	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.59%	Poor
Huang, Wang	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	92.68%	Poor
Zhang, Cai	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.64%	Poor
Guo	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	90.37%	Poor
Liu, Ming	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	90.24%	Poor
Huang, Yang	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	88.89%	Poor
Xu	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	88.68%	Poor
Li	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	82.35%	Poor
Rentsch	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	5.30%	Fair
Hu	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	88.24%	Poor
Wang, Pan	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	87.20%	Poor
Chow (US CDC)	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	No	No	Not applicable	No	96.36%	Poor
Petrilli	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	0.00%	Poor
Dong, Cao	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	88.89%	Poor
Kim	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	82.14%	Poor
Shi, Yu	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	91.79%	Poor
Yang, Yu	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	96.15%	Poor
Argenziano	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	0.00%	Fair
Solis	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	90.62%	Poor
Richardson	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	37.42%	Poor
Fontanet	Yes	Yes	Cannot determine	Yes	No	No	Cannot determine	No	No	No	Yes	No	Not applicable	Yes	0.00%	Poor
Zheng, Gao	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	87.88%	Poor
Liao, Feng	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	99.57%	Poor
Rodriguez	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	0.00%	Poor
Magagnoli	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	85.87%	Poor
Shi, Ren	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	89.55%	Poor
Hadjadj	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	0.00%	Fair
Niedzwiedz	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	0.59%	Fair
Gold (US CDC)	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	94.75%	Poor
Yu, Cai	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	91.58%	Poor
Zheng, Xiong	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	0.00%	Poor
Miyars	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	1.88%	Fair
de la Rica	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	79.17%	Poor
Yin, Yang	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	83.02%	Poor
Gaibazzi	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	0.00%	Fair
Shi, Zuo	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	69.79%	Poor
Cho	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	0.00%	Fair
Allenbach	Yes	Yes	Cannot determine	Yes	No	No	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.42%	Poor
Robilotti	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Yes	No	1.65%	Fair
The Opensafely Collaborative	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	4.16%	Fair
Borobia	Yes	Yes	Cannot determine	Cannot determine	No	No	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	92.95%	Poor
Giacomelli	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	0.00%	Poor
Shah	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	No	Cannot determine	Yes	No	Yes	No	23.73%	Poor
Bello-Chavolla	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	90.06%	Poor
Kolin	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	0.81%	Fair
Lubetzky	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	77.78%	Poor
Goyal	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	94.91%	Poor
Feng	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	90.76%	Poor
Yao	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	96.30%	Poor
Sami	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	0.00%	Poor
Almazeedi	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	0.00%	Poor
Carillo, Vega	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	0.00%	Poor
Yanover	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	0.00%	Fair
Hamer	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	No	No	Not applicable	Yes	0.00%	Fair
Regina	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	95.50%	Poor
de Lusignan	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	Yes	13.44%	Fair
Targher	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	91.74%	Poor
Valenti	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	74.14%	Poor
Feuth	Yes	Yes	Cannot determine	Yes	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	0.00%	Fair
Heng	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	86.27%	Poor
Parrotta	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	2.63%	Fair
Shekar	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	52.00%	Poor
Megia-Vilet	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	93.01%	Poor
Chen, Jiang	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	90.37%	Poor
Jiong	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	94.35%	Poor
Rimland	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	81.82%	Poor
Palaiodimos	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	0.00%	Poor
Ip	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	14.61%	Fair
Heili-Frades	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	11.16%	Poor
Vaquero-Roncer	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	93.15%	Poor
Kim, Garg	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	0.08%	Poor
Wu	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	66.67%	Poor
Hopkinson	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	89.04%	Poor
Shi, Zhao	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	95.05%	Poor
Kimmig	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	63.33%	Poor
Al-Hindawi	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	Cannot determine	Yes	No	Not applicable	No	0.00%	Fair
Basse	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	82.27%	Poor
Freites	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	96.75%	Poor
Alshami	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	82.03%	Poor
Russell	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	No	Cannot determine	Yes	No	Not applicable	Yes	29.25%	Poor
Berumen	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	Yes	0.00%	Poor
Gianfrancesco	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	13.67%	Poor
Li, Long	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Cannot determine	No	94.48%	Poor
Batty	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	88.77%	Poor
Israel	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	Yes	No	Yes	No	Not applicable	Yes	0.00%	Fair
del Valle	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	Yes	No	Cannot determine	Yes	No	Not applicable	Yes	71.16%	Poor
Zuo, Zuo	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	72.73%	Poor
Chaudhry	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Yes	No	85.00%	Poor
Louis	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	54.55%	Poor
Soto-Mota	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	88.00%	Poor
Patel	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	9.62%	Poor
Garibaldi	Yes	Yes	Cannot determine	Cannot determine	No	Cannot determine	Cannot determine	No	No	Cannot determine	Yes	No	Not applicable	No	71.88%	Poor
Docherty	Yes	Yes	Cannot determine	Cannot determine	No	Yes	Cannot determine	Yes	No	Cannot determine	Yes	No	Yes	Yes	29.55%	Poor
Boulware	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	96.71%	NA

a <- table_6_word %>%
  filter(`Overall rating (poor, fair, good)` == "Fair")

nrow(a)

## [1] 19

nrow(table_6_word)-nrow(a)

## [1] 84